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CHAPTER  1 


INTRODUCTION 


1.1  Motivation 

The  study  of  atmospheric  turbulence  is  essential  for  the  understanding  of  many  envi¬ 
ronmental  processes,  such  as  transport  of  heat,  momentum  and  matter.  It  has  major  impact 
on  weather  and  climate  modeling  and  needs  to  be  taken  into  account  for  the  design  and 
operation  of  a  wide  variety  of  man-made  structures,  including  buildings,  bridges,  cars  and 
aircrafts.  Hence  the  measurement  of  turbulence  is  of  high  importance. 

Wyngaard  (1981)  describes  advantages  and  disadvantages  of  several  in-situ  measurement 
devices.  The  described  measurement  devices  include  cup-,  propeller-  and  ultrasonic  anemome¬ 
ter  (sonic).  In  Muschinski  and  Lenschow  (2001),  the  sonic  anemometer  was  highlighted  as 
a  highly  reliable  device,  which  has  had  significant  impact  on  the  research  of  atmospheric 
turbulence.  However,  it  has  to  be  installed  carefully  to  avoid  resonances  that  may  corrupt 
the  measurement  as  pointed  out  by  Siebert  and  Muschinski  (2001). 

Two  sonics  described  by  R.M.  Young  Company  (1999)  were  purchased,  installed  on  tripods 
and  used  to  measure  atmospheric  turbulence  data  in  the  surface  layer.  These  sonics  will  be 
part  of  an  integrated  measurement  system.  This  integrated  system  will  allow  the  user  to 
take  advantage  of  the  capabilities  of  diverse  turbulent  measurement  systems  and  provide 
extended  opportunities  to  study  atmospheric  turbulence. 


1 


1.2  History  and  Current  State 


Early  work  on  turbulence  was  done  by  Taylor  (1938)  and  Simmons  et  al.  (1938).  Kol¬ 
mogorov  (1941a)  and  Kolmogorov  (1941b)  describe  the  energy  dissipation  of  turbulence 
and  the  local  structure  of  turbulence,  respectively.  At  the  same  time,  Obukhov  (1941)  dis¬ 
cussed  the  spectral  energy  distribution  of  a  turbulent  flow. 

Later  the  classical  description  of  surface  layer  turbulence  was  published  by  Monin  and 
Obukhov  (1954),  providing  a  generalization  of  the  law  of  the  wall.  Monin  (1970)  gives 
a  review  on  the  research  of  the  ABL,  including  a  summary  of  what  is  now  known  as  the 
Monin-Obukhov  similarity  theory.  An  experimental  investigation  of  the  profile-functions 
resulting  from  the  Monin-Obukhov  theory  and  a  generalization  of  the  turbulence  spectra 
was  given  by  Kaimal  et  al.  (1972)  based  on  the  so  called  Kansas  1968  experiment. 

The  theory  wave  propagation  through  locally  homogeneous  and  isotropic  turbulence  was 
first  described  by  Tatarskii  (1961)  and  Tatarskii  (1971).  This  was  an  especially  important 
contribution  for  the  later  development  and  understanding  of  atmospheric  and  oceanic  re¬ 
mote  sensing. 

Recently,  the  state  of  atmospheric  turbulence  research  has  been  reviewed  by  Muschinski 
and  Lenschow  (2001).  They  confirmed  the  major  progress  made  in  the  1970s  in  the  field, 
but  also  pointed  out  that  since  then  the  growth  of  the  fundamental  research  has  significantly 
slowed  down.  During  the  workshop  that  was  the  basis  of  their  paper,  they  identified  five 
areas,  that  needed  further  research.  One  of  these  was  In-Situ  Measurements.  The  potential 
of  this  area  was  illustrated  by  the  example  of  the  ultrasonic-anemometer  (sonic)  and  its  im¬ 
pact  on  the  research.  Another  area  was  remote  sensing,  illustrated  by  the  impact  of  Radar 
and  Lidar  Remote  Sensing. 

In  this  context,  it  should  be  noted,  that  while  remote  sensing  techniques  can  probe  a  large 
volume  or  area  of  the  atmosphere,  they  often  suffer  from  a  coarse  resolution  compared  to 
in-situ  techniques.  This  is  because  they  measure  spatial  averages  (path  averages,  area  av¬ 
erages,  or  volume  averages)  rather  than  point  values.  However,  to  calibrate  remote  sensors 
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and  to  confirm  their  proper  operation,  in-situ  measurements  ean  serve  as  a  highly  aceurate 
reference.  Thus  the  eombination  of  both  techniques  is  mutually  beneficial. 

In  addition,  the  emerging  computer,  network,  timing  and  data  acquisition  technology  has 
opened  up  a  multitude  of  fundamental  research  opportunities  on  turbulenee  by  allowing  to 
better  merge  data  and  measurement  systems. 

1.3  Summary  of  Chapters 

This  chapter  was  an  outline  of  atmospheric  turbulence  and  its  importanee.  It  also  gave 
a  summary  of  the  historical  milestones  and  the  eurrent  state  of  turbulenee  research. 
Chapter  two  explains  the  physical  structure  of  turbulenee  based  on  the  easeading  process 
in  section  2.1.  The  following  seetion  2.2  reviews  the  classieal  theory  of  locally  homoge¬ 
neous  and  isotropic  turbulence.  This  includes  an  explanation  of  the  autoeovariance  function 
(ACF),  the  structure  function,  the  power  speetmm  and  the  wavenumber  speetrum.  The  in¬ 
tegral  time  and  the  length  seales  are  diseussed  as  well,  as  the  expeeted  5/3  and  the  2/3  laws 
in  the  inertial  subrange  for  the  power  speetrum  and  the  strueture  funetion,  respeetively.  De- 
seription  of  the  4/3  ratio  between  longitudinal  and  the  transverse  eomponents  of  isotropie 
turbulence  eoncludes  this  seetion.  Seetion  2.3  diseusses  the  Monin- Obukhov  theory  and 
gives  definitions  for  the  atmospherie  boundary  layer  and  the  atmospheric  surface  layer. 
The  impact  of  the  roughness  length  is  diseussed.  As  a  next  step,  the  turbulent  fluxes  and 
the  turbulent  sealing  parameter  are  defined.  Seetion  2.3  ends  with  an  explanation  of  the 
density  stratification  and  its  impact  on  stability  conditions,  the  logarithmie  law  of  the  wall 
and  its  generalization  with  the  eorrection  funetion  for  different  stability  eonditions. 

The  third  chapter  eontains  a  deseription  of  the  measurement  and  data  storage  system.  This 
ineludes  the  measurement  principle  of  the  sonic  in  section  3.1,  providing  a  basic  under¬ 
standing  of  how  it  extracts  temperature  and  velocity  from  the  environment.  Then  the  data 
logging  software  is  described  by  means  of  a  program  strueture  diagram  in  seetion  3.2.  Fi- 
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nally  the  time  stamps  are  validated  in  section  3.3  and  the  file  format  of  the  storage  data  is 
shown  in  section  3.4,  respectively. 

In  the  fourth  chapter,  the  experiment  location  and  the  setup  are  described  in  section  4. 1 .  A 
note  on  the  weather  conditions  during  the  experiment  is  included  in  section  4.2  and  the  co¬ 
ordinate  systems  are  defined  in  section  4.3.  Section  4.4  explains  how  the  data  was  prepared 
for  the  analysis.  The  data  preparation  includes  a  description  of  the  raw  data,  a  discussion 
on  the  averaging  interval,  the  concept  of  data  rotation  between  coordinate  systems  and  data 
detrending.  Finally  an  example  is  given  for  raw,  rotated  and  detrended  data. 

In  section  4.5  the  data  is  analyzed  by  determining  the  noise  floor  and  discussing  the  mea¬ 
surement  accuracy.  Followed  by  a  discussion  of  the  autocovariance  function  and  the  struc¬ 
ture  function  in  time  domain  representation,  the  periodogram  and  the  power  spectrum  are 
used  for  a  spectral  analysis.  Based  on  the  time  and  spectral  domain  description,  the  inertial 
subranges  are  shown  and  the  velocity  and  temperature  structure  parameters  are  estimated. 
Finally  the  length  scale  and  the  time  scale  of  turbulence  are  computed  and  isotropy  is  dis¬ 
cussed. 

The  last  chapter,  five,  briefly  summarizes  the  results  with  view  on  the  theoretical  predic¬ 
tions,  concludes  the  validity  of  the  sonic  measurements  and  gives  a  few  suggestions  for 
further  improvement  for  future  turbulence  measurements  with  the  the  presented  system. 
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CHAPTER  2 


SURFACE  LAYER  TURBULENCE  THEORY 


In  this  chapter,  the  basic  theory  for  the  analysis  of  turbulence  measurements  in  the 
atmospheric  surface  layer  is  summarized.  In  section  2.1  a  qualitative  description  of  the 
structure  of  turbulence  is  stated,  to  provide  an  understanding  for  the  two  components  of  the 
surface-layer  turbulence  theory.  These  are  the  classical  theory  of  locally  homogeneous  and 
isotropic  turbulence,  reviewed  in  section  2.2,  and  the  Monin-Obukhov  theory,  discussed  in 
section  2.3. 


2.1  The  Structure  of  Atmospheric  Ihrbulence 

One  simple  way  to  explain  the  physical  process  of  atmospheric  turbulence  is  to  study 
its  spectrum.  Figure  2.1  shows  schematically  a  typical  turbulence  spectrum.  The  spatial 
wavenumber  k  is  a  (reciprocal)  measure  of  the  size  I  of  the  turbulent  eddies, 

(2.1) 

and  the  wavenumber  spectrum  is  the  energy  density  depending  on  eddy  size.  Through 
Taylor’s  Frozen  Turbulence  Hypothesis  (Taylor,  1938),  which  assumes  that  the  mean  wind 
speed  U  is  significantly  larger  than  the  rotation  velocity  of  the  eddies,  the  equation  relating 
spatial  wavenumber  n  and  frequency  u  can  be  written  as 


u 


(2.2) 
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Figure  2.1.  Spectrum  of  atmospheric  turbulence;  from  Garratt  (1992)  [p.l7] 

In  figure  2. 1  three  different  ranges  are  shown,  the  energy  input  range,  the  inertial  subrange 
and  the  viscous  dissipation  range. 

In  the  energy  input  range  turbulence  is  created  through  surface  stress  resulting  in  overturn¬ 
ing  air  or  through  heating  or  cooling  that  create  convection.  In  this  range  the  wavenumber 
or  frequency  is  low,  which  means  that  the  eddies  are  large.  In  the  atmospheric  surface  layer, 
the  boundary  of  the  energy  input  range  to  the  inertial  subrange  is  typically  on  the  order  of 
meters  to  tens  of  meters.  Eddies  in  the  energy  input  range  break  up  and  enter  the  inertial 
subrange,  where  they  get  smaller  and  smaller  through  a  cascading  process.  During  this 
process,  the  kinetic  energy  is  transferred  down  to  smaller  scales,  but  not  dissipated.  This 
goes  on  until  the  third  region,  the  viscous  dissipation  range,  is  reached.  Here  the  eddies 
get  so  small  that  they  are  converted  to  heat.  The  boundary  between  the  viscous  dissipation 
range  to  the  inertial  subrange  is  typically  on  the  order  of  mm’s.  In  the  atmospheric  surface 
layer,  the  lower  and  upper  boundaries  of  the  inertial  subrange  are  usually  estimated  as 

(2.3) 

2:  Vk 
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where  2:  is  the  height  above  ground  level(AGL)  and 


is  known  as  the  Kolmogorov  scale,  where  v  is  the  kinematic  viscosity  of  air  and  e  is  the 
dissipation  rate  of  turbulent  kinetic  energy  per  mass  unit,  or  simply  the  energy  dissipation 
rate. 

2.2  Theory  of  Locally  Homogeneous  and  Isotropic  Ihrbulence 

Subsection  2.2.1  introduces  the  autocovariance  functions  (ACF)  Bf  (r)  and  Bf  {r),  the 
structure  functions  Dj  (r)  and  Df  (r)  for  random  time  functions  /  (t)  and  one-dimensional 
spatial  fields  /  (x),  respectively.  Here  r  is  the  time  lag  and  r  is  the  spatial  lag  in  the  x- 
direction.  In  Subsection  2.2.2,  the  frequency  spectrum  Wf  (u)  is  expressed  in  terms  of 
Bf  (r)  and  Dj  (r),  respectively,  by  means  of  the  Wiener-Khintchine  relationships.  Cor¬ 
respondingly,  the  one-dimensional  wave-number  spectrum  Vf  (k)  is  expressed  in  terms  of 
Bf  (r)  and  Df  (r)  in  subsection  2.2.3.  The  integral  time  scale  r*,  the  integral  length  scale 
Tj  and  the  spectral  widths  of  the  power  spectrum  Ui  and  the  wavenumber  spectrum  Ki  are 
defined  in  subsection  2.2.4.  Subsection  2.2.5  introduces  the  inertial-range  laws  for  Df  (r) 
and,  by  means  of  Taylor’s  hypothesis,  for  Df  (r).  Subsection  2.2.6  gives  the  inertial-range 
relationships  between  the  structure  functions  and  the  power  spectra  of  the  stream-wise 
velocity  u,  the  lateral  velocity  v,  and  the  vertical  velocity  w  in  the  case  of  locally  homoge¬ 
neous  and  isotropic  inertial-range  turbulence. 

2.2.1  Autocovariance  Function  and  Structure  Function 

Following  Tatarskii  (1971)[p.5],  the  temporal  ACF  Bf{ti,t2)  of  a  random  time  function 
f{t)  is  defined  as 
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(2.5) 


where  ()  denotes  the  ensemble  average  and  [  ]*  denotes  the  conjugate  complex  of  the 
quantity  in  the  brackets.  ti  and  t2  are  two  consecutive  measurement  times.  If  f{t)  is  real 
and  stationary,  then  Bf{ti,t2)  can  be  rewritten  as  a  function  of  the  time  lag  r  only: 

Bf{r)  =  {[f{t  +  r)  -  (/)][/(t)  -  (/)]),  (2.6) 

where 

(/)  =  ifiti))  =  (2.7) 

with 

Bj{t)  =  Bfi-r)  (2.8) 

and 

B/(t)  <  B;(0)  =  a}{t)  (2.9) 

where  aj  is  the  variance  (ACF  at  zero  lag).  That  is,  Bfij)  is  an  even  function  and  attains 
maximum  value  at  r  =  0. 

Stationarity  means  that  the  statistics  of  the  random  process  are  time-invariant.  This  implies 
that  variance  and  mean  of  the  process  are  assumed  to  be  independent  of  time: 

(Tj(t)  =5/(0)  =  const.  (2.10) 

and 

{f  if))  =  {f)  =  const.  (2.11) 

This,  however,  is  not  generally  the  case  for  surface  layer  turbulence,  because  all  statistics 
typically  undergo  a  diurnal  cycle.  But  if  the  observed  time  interval  is  short  enough  com- 
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pared  to  the  diurnal  cycle,  the  random  process  describing  the  surface  layer  turbulence  can 
be  approximated  to  have  a  stationary  first  increment : 

Fy(()  ={/(«  + T)  -  {/(«)  (2.12) 

with 

?/(«)  =  /(*)-(/(«)>.  (2-13) 

which  implies  that  the  ensemble  average  is  assumed  to  vary  only  linearly  with  time.  From 
here,  there  are  two  possible  ways  to  proceed.  The  first  is  to  detrend  the  data,  which  makes 
the  data  stationary  again,  and  to  consider  the  ACF  of  the  detrended  process.  The  second 
approach  is  to  introduce  the  structure  function  D{t  1,12)  as  first  suggested  by  Kolmogorov 
(1941a): 

DAhM)  =  ({[/(ti)  -  (/(ti))]  -  [/(t2)  -  {f{t2))]f)  (2.14) 

For  a  linearly  varying  mean  (as  assumed  with  stationary  increments)  the  structure  function 
becomes  a  function  of  the  time  lag  t  =  ti  —  t2  only: 

Df{ti,t2)  =  Df{ti  -  t2)  =  Df{T)  =  {{[f{t  +  r)  -  {f{t  +  r))]  -  [f{t)  - 

(2.15) 

For  the  more  restrictive  case  of  a  constant  mean 


(/(ti))  =  (/(f2))  =  const.  (2.16) 

we  have 

C/W  =  ({|/(«  +  t) -/(()]}").  (2.17) 

For  an  even  more  restrictive  condition  of  stationarity,  a  relation  between  the  ACF  and  the 
stmcture  function  can  be  written  as: 

£>,(r)  =  2lBj(0)  -  %(t)|.  (2.18) 
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Similarly,  the  spatial  ACF  -B/(r)  and  the  spatial  structure  function  Bf{r)  can  be  expressed: 


Bf{r)  =  {[f{x  +  r)  -  {f)][f{x)  -  (/)])  (2.19) 

and 

D/(r)  =  2|B,(0)  -  B/(r)|,  (2.20) 

where  r  is  the  spatial  lag  in  the  x-direction. 


2.2.2  Power  Spectrum 

The  relationships  between  the  power  spectrum  W (u)  and  the  temporal  ACF  B(t)  can 
be  written  in  terms  of  the  Fourier  transformation  and  the  inverse  Fourier  transformation  as 


1  /■“ 

Wf(u)  =  —  /  e-^‘^^Bf(T)  dr  (2.21) 

27r 


and 


Bf{T)=  /  e^^^Wf{uj)duj. 


(2.22) 


These  two  equations  are  known  as  the  Wiener-Khinchine  relationships.  Inserting  equation 
2.22  into  equation  2.18  gives  a  relationship  between  D{t)  and  W (cu): 


D,(t)  =  2£ 


[1  —  COs{uT)]Wf{u)  duj. 


(2.23) 


Differentiation  with  respect  to  r  leads  to  expressions  of  Wf{uj)  in  terms  of  and 

i9"(r),  respectively: 


1 

27rcc; 


sin(a;r)i9j(r)  dr 


(2.24) 


and 


WM  = 


cos{u}t)Dj{t)  dr. 


(2.25) 
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2.2.3  Wavenumber  Spectrum 


Similar  to  the  preceding  section,  the  Wiener-Khinchine  relationships  between  the  one¬ 
dimensional  wave-number  spectrum  along  the  x-direction  V/(k)  and  the  spatial  ACF  Bf{r) 
can  be  written: 

Vf{K)  =  —  [  e-^^^BAr)  dr  (2.26) 

27r 


and 


Bf{r)=  A^Wf{K)dK. 


(2.27) 


Inserting  equation  2.27  into  equation  2.20  gives  a  relationship  between  Df{r)  and  V/(/t): 


Df(r) 


—  cos{K,r)]Vf{K,)  dn. 


(2.28) 


Differentiation  with  respect  to  r  leads  to  expressions  of  Vf^n)  in  terms  of  and  D'j-{r), 


respectively: 


and 


1  /■“ 

Vf{K)  = -  /  sm{Kr)D'f{r)  dr 

271K  Jq 


1 

271k‘^ 


cos{K,r)D"{r) 


dr. 


(2.29) 


(2.30) 


2.2.4  Integral  Scales  and  Spectral  Widths 

In  Tatarskii  (1971)  [p.  10],  the  integral  time  scale  is  defined  as 


and  the  spectral  width  of  the  power  spectrum  is  defined  as 
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with  the  following  relationship  between  the  two: 


TiUJi  =  271. 


(2.33) 


For  the  integral  length  scale  by  means  of  the  Taylor’s  Hypothesis  we  can  write: 


1  /■“ 


(2.34) 


and  for  the  spectral  width  of  the  wavenumber  spectrum 


1  /‘°° 


(2.35) 


with  the  following  relationship: 


KiTj  =  277. 


(2.36) 


2.2.5  Inertial  Subrange  and  Structure  Parameter 


According  to  Kolmogorov  (1941b)  and  Obukhov  (1941)  Df  (r)  is  proportional  to 


in  the  inertial  subrange 


Df  (r)  = 


(2.37) 


where  Cj  is  referred  to  as  the  structure  parameter  of  /.  From  the  relationship  between 
Vf  (k)  and  D'f  (r)  [Equation  2.29],  we  find 


V,(k)  = 


siiihir  -Cir^~  dr 
3  ^ 


=  — —  /  sin(/T;r)r3  ^dr 

SttkJo  ^  ^ 

-  S--r(2/3)sm(|) 
=  0.1244C'jK-i 


(2.38) 
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Taylor’s  hypothesis  leads  to  a  simple  relationship  between  time  lag  and  stream-wise  spatial 
lag: 

r  =  Ut,  (2.39) 

where  U  is  the  magnitude  of  the  ensemble-averaged  wind  vector.  Inserting  into  (2.37)  gives 

=  (2.40) 


From  the  relationship  between  Wf  (u)  and  D'r  (r)  [Equation  2.24],  we  find 


U'Xw)  = 


27rct; 


sm  UT 


3  •' 


-1 


dr 


re. 

=  — -  /  sm(a;r)r3  dr 

Sttu  Jq 

Cp‘^1^  2  tt 

=  0.1244C'Jf/2/^cu“i. 


(2.41) 


It  should  be  pointed  out  that  equations  2.38  and  2.41  are  integrals  from  zero  to  infinity 
over  the  two-sided  spectrum  Vfp)  and  Wf{u),  respectively.  Because  VjP  and  Wf{u) 
are  even  functions,  they  can  be  rewritten  as  one-sided  spectrum  FfP  and  Sf{u)  in  terms 
of  the  two-sided  spectrum,  respectively: 

Ff{K)  =  2Vf{K)  (2.42) 

and 

Sf{u)  =  2Wf{u),  (2.43) 

where  FjP  and  Sf{u)  are  zero  for  negative  k  and  u,  respectively.  Then  equations  2.38 
and  2.41  can  be  written  as 

Ffp)  =  0.2488C'jK-t  (2.44) 
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and 

Sf{u)  =  0.2mcp^/^u-l ,  (2.45) 

where  equation  [2.44]  agrees  with  Muschinski  et  al.  (2001)  equation  [9],  if  it  is  approxi¬ 
mated  as 

Ff{K)  ^0.249CjK-l.  (2.46) 


2.2.6  Locally  Isotropic  Velocity  Fields 

For  locally  homogeneous  and  isotropic  turbulent  velocity  fields,  Kolmogorov  (1941a) 
predicts  a  4/3-ratio  between  the  directional  and  the  normal  components  in  the  inertial  sub¬ 
range,  as  described  by  equation  [21]  of  Kolmogorov  (1941a): 


[1  -  Rnnjr)]  _  4 

[1  -  Rdd{r)]  3 


(2.47) 


where 


Rrt.n.  (r  ) 


B 

n.n.  (r)  B 

n.n  (r) 


BnniO)  a. 


2 

nn 


(2.48) 


is  the  correlation  function  of  the  velocity  component  in  a  direction  normal  (transverse)  to 
the  r-direction  and 


Rdd(r)  = 


Bdd(r)  Bddir 


dd  I 


-Bdd(O) 


O': 


(2.49) 


dd 


is  the  correlation  function  of  the  directional  (or  longitudinal)  velocity  component.  Noting 
that  the  variance  of  the  transverse  component  and  that  of  the  longitudinal  component, 
must  be  equal  for  isotropy, 

(2.50) 


2  _  2 
Onn  —  Odd 


and  that  the  correlation  functions  Rnn  and  Rdd  are  by  definition  1  at  a  spatial  lag  r  =  0, 


-Rrin(O)  —  RddiS^)  —  1, 


(2.51) 
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equation  [2.47]  can  be  rewritten  as: 


al^[Rnn{^)  -  Rnnir)]  _  4 
'^ddi^ddiO)  —  Rddix)^  3 


With  equation  [2.20]  we  find 


Dr, 


Dddir) 


4 

3’ 


(2.53) 


Inserting  equation  [2.37]  into  equation  [2.53]  a  4/3  ratio  is  also  found  for  the  structure 


parameters: 

DUr)  ^  ^  ^  ^  ^  _  4 

0..(r)  Cl,  Cl  Cl  3' 

where  is  the  structure  parameter  of  the  transversal  vector  component  and  is  the 

structure  parameter  of  the  longitudinal  vector  component.  Here,  is  the  structure  param¬ 
eter  of  the  velocity  vector  component  in  the  x-direction,  and  C/  and  are  the  structure 
parameters  of  the  velocity  components  v  and  w,  where  v  is  the  y-component  of  the  velocity 
vector,  and  w  is  the  2;-component.  In  the  following  we  assume  the  2;-axis  to  point  upward. 


2.3  Monin-Obukhov  Theory 

Subsection  2.3.1  gives  a  definition  of  the  atmospheric  boundary  layer  and  the  surface 
layer.  In  subsection  2.3.2,  the  roughness  length  zq  is  defined.  Equations  for  the  turbu¬ 
lence  fluxes  (latent  heat,  sensible  heat,  momentum  and  temperature)  are  given  in  subsec¬ 
tion  2.3.3.  The  turbulence  scaling  parameters,  describing  the  intensity  of  turbulence,  are 
defined  in  section  2.3.4.  Section  2.3.5  introduces  the  stability  conditions.  Section  2.3.6 
states  the  logarithmic  law  of  the  wall.  Profile  correction  functions,  as  the  central  result  of 
the  Monin-Obukhov  theory  (generalization  of  the  logarithmic  law  of  the  wall)  are  discussed 
in  subsection  2.3.7. 


15 


2.3.1  Atmospheric  Boundary  and  Surface  Layer 


Garratt  (1992)[p.l]  defines  the  atmospherie  boundary  layer  (ABL)  as  “the  layer  of  air 
direetly  above  the  Earth’s  surface  in  which  the  effects  of  the  surface  (friction,  heating, 
cooling)  are  felt  directly  on  time  scales  less  then  a  day  and  in  which  significant  fluxes  of 
momentum,  heat  or  matter  are  carried  by  turbulent  motions  on  a  scale  of  the  order  of  the 
depth  of  the  boundary  layer  or  less.” 

The  surface  layer  of  the  atmosphere,  the  lowest  part  of  the  ABL,  is  defined  as  the  region  in 
which  turbulent  fluxes  may  be  assumed  to  be  constant  as  a  function  of  height  and  the  effect 
of  the  Coriolis  force  due  to  the  rotation  of  the  earth  can  be  neglected.  The  height  of  the 
surface  layer  is  commonly  defined  as  10  %  of  the  ABL  height;  see  figure  2.2. 

2.3.2  Roughness  length 

Another  important  parameter  shown  in  figure  2.2  is  the  roughness  length  zq.  The  rough¬ 
ness  length  characterizes  the  surface  stmcture  and  can  have  values  from  a  fraction  of  a 
millimeters  (smooth  surface)  to  several  meters  (rough  surface,  trees,  urban  area).  A  quan¬ 
titative  definition  of  zq  follows  from  the  logarithmic  law  of  the  wall,  see  subsection  2.3.6. 


2.3.3  Turbulence  Fluxes 

Since  the  momentum-flux,  sensible  heat-flux  and  latent  heat  flux  are  approximately 
constant  in  the  surface  layer  (independent  of  height),  the  surface  layer  is  also  referred  to  as 
the  constant  flux  layer.  The  momentum  flux  (also  Reynolds  stress)  is  given  as 

/ _ 2  _ 2 

Tr  =  PairV  U'W'  +  v'w'  ,  (2.55) 
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The  atmospheric  boundary  layer 


Outer  (Ekman)  layer 
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Figure  2.2.  Structure  of  atmospheric  boundary  layer;  from  Garratt  (1992)[p.2] 
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where  u'w'  is  the  covariance  of  the  horizontal  velocity  u'  with  the  vertical  velocity  w'  and 
where  v'w'  is  the  covariance  of  the  horizontal  lateral  velocity  v'  with  w' .  pair  is  the  density 
of  air.  The  sensible  heat  flux  Hs  is  given  by 

Hs  PairC^p^  ^  )  (2.56) 

where  Cp  is  the  specific  heat  capacity  of  air  at  constant  pressure  and  w'9'  is  the  covariance 
of  w'  with  the  temperature  6'.  The  latent  heat  flux  is  defined  as 

Hl  =  \wU!'q’,  (2.57) 

where  is  the  specific  heat  of  vaporization  of  water,  and  w'q'  the  covariance  of  w'  with 
the  specific  humidity  q'.  In  addition,  a  temperature  flux  Q  can  be  denoted  as 

Q  =  (2.58) 

Pair  Cp 

which  is  sometimes  referred  to  as  the  kinematic  heat  flux. 


2.3.4  Turbulence  Scaling  Parameters 

The  following  four  scaling  parameters  are  defined  for  the  surface  layer.  These  are: 
velocity  scale  (also  known  as  the  friction  velocity) 

=  (W%^^)^/^  (2.59) 


the  temperature  scale 


Q  w'9' 


(2.60) 
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the  Monin- Obukhov  length 


uj6  uj^e  9 

*  “  KgQ  ~  Kg^  ~  Kg  9,  ^ 

and  the  height  2:  AGL.  9  denotes  the  mean  temperature  in  the  surface  layer,  K  the  Karman 
constant  and  g  is  the  constant  of  gravity. 

2.3.5  Stability  Conditions 

Stability  conditions  can  be  divided  into  the  cases  of  unstable,  neutral  and  stable  strat¬ 
ification.  The  dimensionless  height  ^  gives  a  compact  description  of  the  three  cases,  in 
which  gives  negative  values  for  unstable  stratification,  values  approaching  zero  for  neu¬ 
tral  stratification  and  positive  values  for  stable  stratification. 


2.3.6  Logarithmic  Law  of  the  Wall 

Equation  [13]  in  Monin  and  Obukhov  (1954)  states  the  gradient  of  the  mean  wind  speed 
with  respect  to  height  2:  (mean  wind  speed  height  profile)  near  the  surface  as: 


dU  {z) 
dz  Kz 


(2.62) 


Similarly  the  gradient  of  the  mean  temperature  with  respect  to  height  (mean  temperature 
height  profile)  can  be  written  as: 


d9{z)  9^ 


dz  z 

Integration  of  2.62  and  2.63  with  respect  to  height  2:  leads  to 


(2.63) 


U{z) 


(2.64) 
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and 


0(z)  =  9,  In  ,  (2.65) 

which  are  the  logarithmic  law  of  the  wall.  Here  equation  [2.64]  is  identical  with  equation 
[18]  in  Monin  and  Obukhov  (1954).  Equations  [2.64]  and  [2.65]  give  the  mean  values  of 
wind  speed  and  temperature  at  a  given  height  2;. 


2.3.7  Profile  Correction  Functions 

In  Monin  and  Obukhov  (1954)  equations  [29']  and  [30']  state  a  more  general  form  of 
the  equations  [2.62]  and  [2.63]  including  all  three  possible  density  stratifications: 


dz 


Y'z 


(2.66) 


and 


d6{z)  _  9^  ^  f  z 


(2.67) 


dz  z 

where  01  is  the  correction  function  of  the  height  profile  of  the  mean  wind  speed  and  02  is 
the  correction  function  of  the  height  profile  of  the  temperature.  Thus  the  Monin-Obukhov 
theory  is  a  generalization  of  the  logarithmic  law  of  the  wall.  Equations  [2.66]  and  [2.67] 
can  be  integrated  with  respect  to  2:  and  lead  to  the  mean  wind  speed  and  temperature  at  a 
given  height  2:  respectively: 


and 


(2.68) 


9iz)  -  9{zo) 


(2.69) 
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where  j  and  ^3/2  yT~j  correction  functions  for  the  integrated  profiles,  corre¬ 

sponding  to  01  and  02-  From  equation  [3.34]  [p.  53]  in  Garratt  (1992)  it  can  be  seen  that  a 
relationship  between  4/  and  0  can  be  expressed  as: 


For  small  values  of  0  j  an  analytical  form  can  be  written  as: 


(2.70) 


(2.71) 


where  01  is  are  the  coefficient  of  a  power  law  expansion  (neglecting  terms  of  order  higher 
than  one)  and  equation  [2.71]  is  identical  to  equation  [42]  in  Monin  and  Obukhov  (1954). 
Inserting  equation  [2.71]  into  equation  [2.70]  gives: 


*  (t)  -  -ft  t-  (2.72) 

Inserting  equation  [2.72]  into  equation  [2.68]  and  equation  [2.69],  results  in  equations  for 
the  mean  wind  velocity  U{z)  and  the  mean  temperature  9{z)  as  given  by  equations  [43a] 
and  [43b]  in  Monin  and  Obukhov  (1954)  respectively: 


U{z)  = 


ln(-)+0if 
zqJ  L* 


(2.73) 


known  as  the  log-linear  wind  law,  and 


e{z)  -  e{zo)  =  ^ 


In  (-)+02f 
\zoJ 


(2.74) 


where  02  is  another  coefficient.  Figure  2.3  shows  qualitatively  how  the  different  stability 
conditions  resulting  from  the  density  stratification  impact  the  wind  profile. 
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CHAPTER  3 


SYSTEM  DESCRIPTION 


This  chapter  explains  the  system  used  for  collection  and  storage  of  wind  velocity,  tem¬ 
perature  and  time  stamps.  The  measurement  principles  of  the  sonic  anemometer  and  ther¬ 
mometer  are  described  in  section  3.1.  Followed  by  the  explanation  of  the  data  logging 
system  in  section  3.2,  which  receives  the  data  from  the  sonic,  adds  a  time  stamp  and  stores 
it,  section  3.3  provides  an  analysis  of  the  accuracy  of  the  time  stamps.  In  section  3.4  the 
format  of  the  time  stamps,  data  and  the  debug  information  is  shown. 

3.1  Sonic  Anemometer  and  Thermometer  Measurement  Principle 

The  sonic  anemometer  relies  on  a  simple  and  efficient  measurement  method,  which 
will  be  summarized,  based  on  Kaimal  and  Finnigan  (1994)[p.247].  In  the  given  sonic 
anemometer  from  R.M.  Young  Company  (1999),  three  measurement  paths  of  length  d 
crossing  through  the  same  measurement  volume  from  different  directions  are  defined,  as 
shown  in  figure  3.1.  An  ultrasonic  transducer  is  placed  at  both  ends  of  each  of  these  paths. 
This  transducer  can  send  and  receive  an  ultrasonic  sound  wave  of  a  specific  frequency,  as 
seen  in  figure  3.2.  The  time  ti  for  the  ultrasonic  sound  wave  to  travel  from  one  transducer 
to  its  counterpart  can  be  stated  as: 


h  = 


(3.1) 
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Figure  3.1.  Sonic  anemometer  measurement  paths  and  sample  volume  (Pieture:  R.M. 
Young  Company  (1999)) 


- ►  c  +  Vd 

I  C  -  Vd  '  t 

- d - ^ 
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Figure  3.2.  Ultrasonie  transdueer  and  measurement  path  from  Kaimal  and  Finnigan  (1994) 
[p.247] 
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where  c  is  the  speed  of  sound  and  Vd  is  the  wind  velocity  along  the  measurement  path.  Sim¬ 
ilarly  the  time  t2  for  an  ultrasonic  wave  traveling  in  the  opposite  direction  can  be  described 
as: 


d 

c+Vd 

Equations  [3.1]  and  [3.2]  can  be  rewritten  as: 


(3.2) 


c-Ki 


ti 


(3.3) 


and 

c+Vd  =  ^.  (3.4) 

t2 

Subtracting  equation  [3.3]  from  [3.4]  gives  the  wind  velocity  along  the  path  as  Vd'. 


Vd 


d  n 
2 


(3.5) 


This  gives  the  wind  velocity  Vd  along  the  path  d  only  in  terms  of  the  time  ti,  the  time  t2  and 
the  path  length  d.  Note  that  for  the  sound  velocity  the  following  empirical  equation  exists: 


=  403T(1 +  0.32-),  (3.6) 

P 

where  T  is  the  air  temperature  (in  K),  e  the  water  vapor  pressure  and  p  the  air  pressure.  In 
Kaimal  and  Finnigan  (1994)  equation  [3.6]  is  approximated  as: 


=  403T(1  +  0.32-)  ^  403T(1  +  0.38-)  =  403T^,  (3.7) 

p  p 

with  as  the  virtual  temperature.  Taking  the  sum  of  equation  [3.3]  and  [3.4]  results  into: 


2c  _  1  1 

d  ti  t2 


(3.8) 
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Equation  [3.8]  can  be  rewritten  as 


(3.9) 


or 


(f 


Equating  equation  [3.7]  and  equation  [3.10]  gives: 


n  ly 

1612  \ti  ’ 


(3.10) 


(3.11) 


which  is,  again,  only  in  terms  of  the  of  the  time  ti,  the  time  t2  and  the  path  length  d. 
This  shows  that  measuring  the  travel  times  ti  and  t2  is  sufficient  to  compute  the  wind 
vector  components  and  the  virtual  temperature.  In  the  given  sonic  this  is  done  by  a  micro¬ 
controller  integrated  into  the  sonic,  which  measures  the  transition  times,  applies  several 
corrections  (e.g.  wake,  cross  wind,  transducer  shadowing)  and  computes  the  wind  vec¬ 
tor  components  and  the  virtual  temperature.  This  data  is  then  provided  through  a  serial 
asynchronous  RS232  interface. 


3.2  Data  Logging  System 

The  serial  interface  of  the  sonic  described  in  subsection  3.1  was  connected  to  a  serial- 
to-universal  serial  bus  (USB)  converter  for  each  of  the  two  sonics  used.  And  these  were 
attached  to  a  standard  notebook  computer  running  Linux. 

In  Linux  any  hardware  or  software  device  or  data  unit  is  represented  by  a  file.  In  the  case 
of  the  serial  ports  provided  via  the  USB  port  these  files  are  /dev/ttyUSBO  and  /dev/ttyUSBl 
respectively.  In  this  way  the  interaction  with  these  devices  is  provided  in  a  similar  way  as 
interaction  is  provided  with  data  units  (normal  files).  So  for  data  units  the  main  actions  are 
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Allocate  program  variables  and  resources 

Open  file  descriptors  (FD)  devices  to  sonic  1  and  2 

Clear  and  initial  new  FD  settings 

Save  old  serial  port  settings 

Initial  new  serial  port  settings 

Write  new  FD  and  serial  port  settings 

Open  files  for  data  and  time  stamps  storage 

Wait  for  a  dummy  read  from  sonic  1 

Wait  for  a  dummy  read  from  sonic  2 

For  index=0  to  infinity 


Active  Sonic  ? 

Sonic  1 

Sonic  2 

Read  data 

Read  data 

Read  time  stamp 

Read  time  stamp 

Write  time  stamp 

Write  time  stamp 

Write  data 

Write  data 

Clear  sonic  1 

Clear  sonic  2 

Activate  sonic  2 

Activate  sonic  1 

Figure  3.3.  Program  structure  of  data  logging  program 
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to  read  or  write  data  or  to  seek  a  position.  For  devices  another  action  is  possible,  this  action 
is  the  configuration  of  the  device  (also  called  setup  control  or  initialization).  This  action 
usually  specifies  the  way  the  device  operates  (e.g.  how  it  handles  data).  However  once  this 
phase  is  finished,  the  device  can  be  treated  very  similar  to  a  data  unit. 

In  figure  3.3,  the  logical  flow  of  the  data  logging  and  time  stamping  program  is  shown, 
starting  with  allocating  variables  and  other  resources  needed.  The  second  step  is  to  open 
the  serial  port  devices  using  a  file  descriptor  (FD)  list.  A  FD  in  this  concept  is  essentially  a 
program  internal  universal  pointer  to  any  device  file.  This  FD-list  is  essential  for  defining 
a  set  of  devices.  Once  opened,  these  FD’s  are  reset  to  a  defined  state  by  emptying  their 
buffers,  in  step  three.  As  the  next  step  the  old  configuration  of  the  serial  port  is  saved  as  a 
precaution  and  for  step  five  the  new  desired  configuration  of  the  serial  ports  is  assigned  to 
variables.  As  step  six  the  configurations  of  the  FD’s  and  the  serial  ports  are  written,  which 
concludes  the  setup  process. 

Now,  as  the  devices  are  ready  to  be  read,  the  program  does  a  dummy  read  on  each  of  the 
two  devices,  disregards  this  data  received  and  clears  the  buffer.  This  action  ensures  that 
the  program  is  now  synchronized  with  the  sonics  and  does  not  receive  incomplete  data 
messages.  As  of  this  point  the  second  sonic  was  the  last  device  read  and  the  first  sonic 
is  assumed  to  provide  data  next  (to  be  set  active  in  the  FD  list).  Once  the  activity  signal 
is  detected,  the  data  is  read  from  the  device  and  a  time  stamp  is  generated.  Both  are  then 
stored  in  the  data  unit.  The  final  step  before  getting  ready  to  receive  data  again  is  now  to 
clear  the  device  again  and  to  activate  the  second  sonic,  which  is  expected  to  receive  data 
next. 

Then  the  whole  process  is  repeated  for  the  second  sonic  by  receiving  data,  creating  the 
time  stamp,  storing  data  and  time  stamp  and  switches  back  to  the  first  sonic.  This  process, 
known  as  ping-pong  technique,  continues  until  an  interruption  signal  is  given  by  sending 
CTRL-C  to  the  program.  After  this  signal,  a  signal  handler  will  be  called.  This  is  the  only 
way  to  exit  the  infinite  loop.  The  signal  handler  then  closes  the  data  units,  restores  the  old 
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serial  port  configurations  and  terminates  the  program. 

Note  that  this  last  step  is  not  shown  in  figure  3.3,  since  it  is  not  part  of  the  main  program, 
but  merely  a  way  to  terminate  the  data  logger  correctly  at  any  given  time. 

3.3  Time  Stamp  Validation 

Section  3.2  pointed  out  that  not  only  the  data  but  also  time  stamps  associated  with  the 
data  are  stored.  It  is  important  to  determine  the  validity  of  these  time  stamps,  and  whether 
they  are  a  reliable  representation  of  time.  Ideally  the  time  stamps  should  be  equidistant  in 
time,  resulting  in  a  linear  increase  of  recorded  time  as  a  function  of  the  time-stamp  index. 
In  the  case  of  interleaved  timestamps  from  two  sensor  data  streams  (as  it  is  the  case  for  the 
two  sonics),  the  offset  between  the  time  stamps  with  even  and  odd  index  should  be  constant, 
which  implies  that  the  sampling  rates  of  the  two  sensors  are  exactly  identical  (e.g  derived 
from  the  same  reference  oscillator).  This  furthermore  assumes  that  the  time  for  both  sonics 
from  measurement  of  the  data  till  time  stamp  generation  and  storage  is  constant. 

For  the  further  analysis  a  few  terms  have  to  be  defined,  which  are  raw  time  stamp,  cor¬ 
rected  time  stamp  and  linear  time  stamp.  Note  that  the  time  stamp  is  originally  created  as 
described  in  section  3.4  for  every  data  point.  For  the  following  analysis  the  time  stamps 
are  converted  to  seconds  counted  from  the  beginning  of  the  experiment.  This  procedure 
simplifies  the  analysis  significantly. 

The  simplest  time  stamp  is  a  raw  time  stamp,  which  is  the  time  stamp  as  it  is  stored  by 
the  data  logging  system.  In  figure  3.4  in  the  upper  and  center  diagram  the  first  hundred 
raw  time  stamps  for  the  first  and  the  second  sonic  are  plotted  as  solid  lines  against  their 
index.  The  spikes  shown  in  these  time  stamps  are  caused  by  an  operating  system  adjusting 
mechanism,  which  tries  to  compensate  for  the  PC-clock  drift. 

A  corrected  time  stamp  means  that  an  algorithm  has  checked  whether  the  time  between  two 
consecutive  samples  is  an  integer  number  multiple  of  the  sampling  interval  or  not.  In  the 
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2005-10-17  17:29:52-2005-10-17  17:29:55 


Figure  3.4.  Raw  and  corrected  time  stamps  from  sonic  1  and  2  data  and  time  difference 
over  a  time  of  4s. 
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case  of  a  negative  result,  the  erroneous  time  stamp  will  be  corrected  by  replacing  it  with  the 
directly  preceding  time  stamp  plus  the  nearest  integer  number  multiple  of  the  sampling  in¬ 
terval.  Note  that  this  algorithm  allows  for  a  small  correct  range  around  the  integer  number 
multiple  of  the  sampling  interval  for  the  time  between  samples,  rather  than  one  value.  This 
is  needed  beeause  there  is  always  small  jitter  in  the  time  stamping  interval.  Therefore  this 
method  does  not  eorreet  for  any  drift  and  might  eontain  jumps  in  time  because  it  eheeks 
for  the  integer  multiple  of  the  sampling  interval,  rather  than  for  the  sampling  interval.  In 
figure  3.4  the  dashed  line  shows  the  eorrected  time  stamps  corresponding  to  the  raw  time 
stamps  plotted  against  their  index.  It  ean  be  seen  that  spikes  from  the  raw  time  stamps  were 
successfully  removed. 

Linear  time  stamping  means  that  the  time  stamps  are  reeomputed  by  taking  the  last  and 
the  first  time  stamp  from  the  whole  data  set,  and  that  a  linear  estimation  between  them  is 
computed.  Then  all  values  between  the  first  and  the  last  time  stamp  are  replaced.  Note 
that  use  of  this  teehnique  ean  introduce  additional  errors  compared  to  the  true  time  if  the 
original  time  series  eontains  a  drift.  Otherwise  the  resulting  time  stamps  will  be  perfeet,  as 
missing  samples  are  assumed  not  to  exist. 

The  bottom  diagram  of  figure  3.4  shows  the  differenees  between  the  two  sonies  for  the 
raw  and  for  the  eorrected  time  stamps.  It  ean  be  seen  that  the  difference  of  the  eorrected 
time  stamp  is  close  to  zero  and  that  the  difference  between  raw  time  stamps  is  a  small  value 
compared  to  the  sampling  interval.  For  the  latter  one,  however,  this  small  offset  is  inereased 
when  the  elock  adjustment  spike  is  observed.  This  seems  to  be  the  ease  beeause  the  cloek  is 
adjusted  mueh  faster  then  normally.  However,  even  this  increased  offset  is  small  eompared 
to  the  time  between  time  stamps,  so  that  no  large  relative  error  is  generated  between  time 
stamps  of  both  sonies. 

Since  figure  3.4  only  shows  a  very  small  time  period,  not  all  possible  errors  might  have 
occurred.  Figure  3.5  shows  the  relative  error  between  the  sonies  for  all  three  available  time 
stamps  over  the  total  experiment  time.  It  ean  be  seen  that  there  are  many  spikes  in  the 
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Ti-T2[s] 


Time  [hours] 

—  Raw 

—  Corrected 

—  Linear 


Figure  3.5.  Relative  error  between  both  sonics  for  linear,  corrected  and  raw  timestamps. 
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raw  time  stamps,  significantly  reduced  number  of  spikes  in  the  corrected  time  stamps  and 
none  in  the  linear  time  stamps.  Assuming  that  few  spikes  in  the  corrected  time  stamps  are 
caused  by  missing  samples,  it  can  be  justified  to  assume  that  most  of  the  data  shown  has  no 
missing  samples.  Based  under  this  assumption  the  problem  of  ignoring  missing  samples, 
when  using  the  linear  time  stamps,  should  be  negligible. 

Here  two  more  possible  errors  are  explained.  Note  that  all  effects  are  superimposed  on 

2005-10-17  17:29:52  -  2005-10-20  13:14:43 


Figure  3.6.  Drift  of  PC  clock  over  66  h. 

each  other.  One  of  the  additional  errors,  which  can  appear  if  two  sonics  are  not  referenced 
to  the  same  oscillator  (the  sampling  frequencies  might  not  be  exactly  the  same),  is  that  a 
phase  error  can  accumulate.  This  can  switch  the  order  of  sampling  of  the  sonics  with  re¬ 
spect  to  each  other.  Thus  the  data  logging  program  described  in  section  3.2  will  lose  one 
data  point,  since  it  expects  a  value  from  the  other  sonic  to  be  sent.  However  this  error  can 
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be  tracked  by  finding  the  difference  between  two  consecutive  samples,  and  detect,  when 
the  difference  is  twice  the  sampling  interval.  Note  that  this  error  still  assumes  a  perfect 
PC-clock  oscillator. 

The  PC-clock  oscillator  drifts  in  general  due  to  temperature  changes.  In  figure  3.6  such  a 
drift  is  shown  for  the  whole  experiment  time  period.  Figure  3.6  shows  the  difference  be¬ 
tween  the  corrected  time  stamps  and  the  linear  time  stamps.  Interestingly,  a  kind  of  diurnal 
cycle  can  be  seen  in  the  drift,  which  confirms  the  dependency  of  the  time  stamp  drift  with 
temperature.  The  spikes  are  likely  resulting  from  missing  samples,  as  explained  before. 
Note  that  the  drift  of  the  PC-clock  over  the  experiment  time  is  on  the  order  of  a  second. 
When  using  the  two  sonics  alone  with  the  same  computer  as  a  data  logger,  it  can  also  be 
seen  that  the  relative  time  between  the  two  sonics  is  stable.  But  the  drift  might  become  a 
problem  in  reference  to  external  systems. 

3.4  File  Format 

The  data  and  the  time  stamps  are  stored  to  a  data  unit  on  a  mass  storage  device. 
The  data  unit  is  a  plain  text  file  (also  known  as  ASCII  file).  The  two  files  containing 
the  two  data  stream  from  sonic  1  and  2  are  named  SONIC  J  J20051017205809.DAT  and 
SONICJ2J20051017205809.DAT,  respectively.  The  time  stamp  included  in  the  filename  is 
the  start  time,  day,  month  and  year  of  the  data  record. 

The  files  contain,  for  each  text  line  in  this  order,  a  time  stamp,  as  shown  in  table  3.1,  the 
data  as  shown  in  table  3.2  and  debugging  information,  as  shown  in  table  3.3.  Note  that  for 
the  main  data  set  used  in  the  data  analysis  (2005 1017T175952  eastern  daylight  time  (EDT) 
to  2005 1020T1 15952  EDT)  only  the  bold  type  values  were  recorded,  to  reduce  the  file  size. 
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Year 

Month 

Day 

Hour 

Minute 

Seconds 

2005 

07 

20 

14 

31 

15.355233 

Table  3.1.  Time  stamp  format 


u 

V 

w 

Azimuth 

Elevation 

3d-Speed 

Speed  of  Sound 

Temperature 

|H| 

■HI 

[°] 

[°] 

■HH 

[-] 

•-  ,s  -* 

[°C'] 

-0.86 

0.16 

-0.20 

281.0 

-12.9 

0.91 

349.84 

30.57 

Table  3.2.  Data  format 


Error  Flag 

Internal  Voltage 

0 

26.5 

Table  3.3.  Debug  information 
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CHAPTER  4 


EXPERIMENT  SETUP  AND  DATA  ANALYSIS 


4.1  Experiment  Location  and  Setup 

To  test  the  turbulenee  measurement  system  deseribed  in  the  preeeding  ehapter,  an  ex¬ 
periment  was  performed  at  the  Tilson  farm  faeility  of  the  Mierowave  Remote  Sensing  Lab 
in  Amherst,  MA.  Timestamps,  temperatures  and  the  three  wind  veetor  eomponents  were 
reeorded  for  both  of  the  sonies  from  the  17th  of  Oetober  2005  at  17:29:52  eastern  daylight 
time  (EDT)  until  the  20th  of  Oetober  2005  at  13:14:42  EDT. 

Eigures  4. 1  and  4.2  show  the  setup  of  the  sonies  on  top  of  standard  light  stand  tripods.  The 
sonies  were  stabilized  using  nylon  strings  attaehed  to  the  sonie  frame  and  to  tent  pegs  in 
the  ground.  The  measurement  height  was  Zi  =  3.8  m  for  one  sonic  and  Z2  =  3.0  m  for  the 
other.  The  centers  of  the  tripod  bases  were  separated  by  d  =  1.8  m.  It  should  be  noted  that 
this  was  only  a  tentative  setup,  since  the  two  professional  20-foot  meteorological  tripods 
ordered  were  not  available  in  time.  In  the  future,  the  professional  system  will  be  used 
because  it  is  more  robust,  allows  larger  measurement  heights  and  can  facilitate  multiple 
sensors  on  a  single  tower. 

4.2  Weather  Conditions 

Another  factor  influencing  the  data  is  the  weather.  Analysis  of  the  cloud  coverage 
information  from  www.weather.com  for  the  given  time  period  from  the  17th  of  October 
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Figure  4.1.  Experiment  setup  at  Tilson  farm  faeility  (Pieture:  Yonghun  Cheon) 


Figure  4.2.  Close  up  of  sonies  (Pieture:  Yonghun  Cheon) 
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2005  at  17:29:52  EDT  till  the  20th  of  October  2005  at  13:14:42  EDT,  shows  that  it  was 
between  cloudy,  mostly  cloudy  and  sometimes  partly  cloudy.  There  was  no  precipitation. 


4.3  Coordinate  Systems 

The  sonic  can  be  set  to  provide  the  wind  vector  components  u  (horizontal),  v  (lateral) 
and  w  (upward)  in  a  rectangular  coordinate  system  referenced  to  its  frame,  as  shown  in 
figure  4.3.  However  the  most  useful  coordinate  system  for  the  application  of  turbulence 
theory  is  a  stream-wise  one,  since  it  allows  for  simplification  of  the  theory.  In  the  rectan¬ 
gular  stream-wise  coordinate  system  the  unit  vector  of  the  wind  component  u  is  the  unit 
vector  of  the  mean  wind  vector.  The  data  is  converted  from  the  sonic  coordinate  system  to 
the  stream-wise  coordinate  system  by  rotation,  which  is  described  in  subsection  4.4.3. 


W 


* 

Figure  4.3.  Sonic  coordinate  system  (Picture:  R.M.  Young  Company  (1999)) 
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4.4  Data  Preparation 


4.4.1  Raw  Data 

Raw  atmospheric  data  can  be  written  in  the  form  of  a  mean  quantity  /  plus  a  fluctuating 
part  /',  also  known  as  Reynolds  notation: 


/  =  /  +  /'• 


(4.1) 


For  a  wind  vector  decomposed  into  the  rectangular  sonic  coordinate  system  the  resulting 
equations  are  shown  in  table  4.1. 


Wind  vector  component 

Total  Component 

Mean 

horizontal 

u  =  u  +  u' 

horizontal  lateral 

V  =  V  +  v' 

vertical  upward 

w  =  w  +  w' 

Table  4.1.  Wind  vector  components  in  Reynolds  notation 


4.4.2  Averaging  Interval 

To  resolve  mean  quantities  of  velocity  and  temperature,  an  averaging  interval  has  to  be 
chosen.  It  has  to  take  two  important  parameter  into  account.  First  of  which  is  the  variance 
of  the  wind  fluctuation  around  the  mean.  The  second  one  is  the  time  scale  of  the  diurnal 
cycle,  which  is  24  hours,  but  can  have  significant  changes  on  the  order  of  hours.  As  a  part 
of  the  Reynolds  notation  in  table  4.1,  the  mean  value  of  a  stationary  random  time  function 
f{t)  was  estimated  as  time  average  f{t)  instead  of  ensemble  average  (/(f)).  However  if  a 
time  average  is  taken  over  a  long  enough  time,  it  converges  to  the  ensemble  average  (/(f)) : 

1  ^ 

fii)  =  ifii))  =  lim  (4.2) 

N^CO  iV  i  o 

i=l 
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where  is  the  sampling  interval,  i  is  the  sample  number  and  N  the  total  numbers  of 
samples  taken  into  account.  Using  an  infinite  sum  of  a  time  series  however  is  impractical 
and  the  law  of  large  numbers  allows  for  truncating  of  this  sum  if  the  number  of  averaged 
independent  samples  Nind  =  —  is  large.  Here  Tj  is  the  integral  time  scale  defined  in 
equation  [2.31]  and  T  =  NTg  is  the  total  time  taken  into  account  for  the  calculation  or 
simple  averaging  time.  A  sample  is  called  independent  in  this  context  if  it  is  spaced  further 
then  the  integral  time  scale  r*  away  from  the  preceding  sample.  In  other  words,  a  sample 
is  dependent  on  a  preceding  sample,  if  the  spacing  between  them  is  less  then  the  integral 
time  scale  r*. 

In  order  to  obtain  a  large  number  of  independent  samples,  the  averaging  time  T  therefore 
must  be  much  larger  then  the  integral  time  scale  Xj: 

T  »  Ti.  (4.3) 


The  law  of  large  numbers  can  be  written  as: 


2  A  2 


(4.4) 


where  ajr  is  the  variance  of  the  time  average  (with  the  averaging  time  T).  To  get  a  real¬ 
istic  understanding  of  the  range  of  expected  reduction  of  the  variance  of  the  estimation  of 
the  mean  value  of  the  velocity  or  temperature,  real  data  has  to  be  inserted  into  equation 
[4.4].  From  velocity  measurements  shown  later  in  section  4.5  the  variance  for  one  given 
time  series  was  estimated  as  =  1.20(y)^  and  the  integral  time  scale  was  estimated  as 
Ti  ~  5s.  Note  that  the  averaging  time  is  a  value  which  has  to  be  a  trade-off  between  two 
demands.  First  there  is  the  demand  for  a  larger  averaging  time  to  reduce  the  variance  of 
the  estimation  of  mean  values.  On  the  other  hand,  we  assume  stationarity,  which  does  not 
hold  for  longer  averaging  times,  which  are  a  significant  fraction  of  the  diurnal  cycle.  In 
Kaimal  and  Finnigan  (1994)  [p.256]  a  maximum  averaging  time  is  suggested  as  about  as 
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T  Ri  3600s,  used  in  the  following  estimations.  Inserting  these  values  gives  the  standard 
deviation  of  the  stationary  random  time  function  of  the  velocity  component  u  as: 


-^1,20(™P  =  0.041™ 

3600s  ^  s’  s 


(4.5) 


Another  way  to  express  this  relationship  is  by  the  number  of  independent  samples,  which 
can  be  calculated  for  the  given  example  as: 


^ind 


T 


Ti 


3600s 

5s 


720. 


(4.6) 


Table  4.2  gives  results  for  the  number  of  independent  samples  and  accuracies  for  different 
averaging  periods  and  a  fixed  integral  time  scale.  As  seen  here,  the  user  of  the  data  has  to 
decided,  what  is  an  acceptable  au,  which  then  sets  the  required  averaging  time. 


Averaging  Time  [s] 

#  of  Independent  Samples 

au  in  %  of  au 

60 

12 

28.9 

600 

120 

20.4 

3600 

720 

3.7 

Table  4.2.  Number  of  independent  samples  and  accuracy  for  integral  time  scale  r*  =  5s 


4.4.3  Data  Rotation 

In  section  4.3  the  sonic-  and  stream-wise  coordinate  systems  were  introduced.  One 
commonly  used  method  to  convert  data  from  the  sonic-  to  the  stream-wise  coordinate  sys¬ 
tem  is  the  three-angle-rotation,  which  is  described  in  Wilczak  et  al.  (2001)  and  summarized 
here.  Figure  4.4  shows  a  visualization  of  the  rotations  involved.  The  first  step  is  to  deter- 
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mine  the  angle  a  from  the  mean  wind  eomponent  u  in  the  x-direetion  and  the  mean  wind 
eomponent  w  in  the  2;-direetion: 

a  =  arctan  .  (4.7) 

The  angle  a  ean  be  used  to  to  eompute  the  rotation  matrix  B: 

cos(q;)  0  sin(Q;) 

B=  0  10,  (4.8) 

—  sin(Q;)  0  cos(q;) 

which  converts  the  raw  wind  vector  component  u  into  the  intermediate  wind  vector  com¬ 
ponent  uj  in  the  intermediate  X/-direction  and  the  raw  wind  vector  component  w  into  the 
intermediate  wind  vector  component  wi  in  the  intermediate  ^/-direction: 

u 

B  V  (4.9) 

w 

Note  that  the  raw  wind  vector  component  v  was  not  effected  by  this  step.  This  operation 
has  now  corrected  the  pitch,  the  next  step  is  to  compute  the  angle  /9  to  correct  the  roll  angle: 

(3  =  arctan  ^ ,  (4.10) 

where  wi  is  the  mean  of  the  intermediate  wind  vector  component  wi.  Following  the  rotation 
matrix  C  can  be  computed: 

1  0  0 

C=  0  cos(/9)  —  sin(/9)  •  (4.11) 

0  sin(Q;)  cos(/9) 
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and  applied  similar  to  the  rotation  matrix  B: 


" 

~ 

Ui 

Ui 

Vi 

=  C 

V 

w' 

Wi 

(4.12) 


which  left  uj  unchanged  and  converted  the  raw  wind  vector  component  v  into  the  inter¬ 
mediate  wind  vector  component  vj  in  the  Ij-direction  and  the  intermediate  wind  vector 
component  wj  into  the  final  wind  vector  component  w'  in  the  Z'-direction  in  the  stream- 
wise  coordinate  system.  Finally  the  yaw  angle  7  can  be  computed  from  the  mean  of  the 
intermediate  wind  vector  component  vj  and  from  the  mean  of  the  intermediate  wind  vector 
component  by  using 

7  =  arctan  •  (4.13) 

The  third  rotation  matrix  can  then  be  found  by  using 


D 


003(7)  —  sin(7)  0 
sin(7)  003(7)  0 

0  0  1 


(4.14) 


Applying  matrix  D  then  gives  the  two  missing  final  wind  vector  components  u'  in  the  X'- 
direction  and  v'  in  the  y'-direction  respectively: 


~ 

~ 

u' 

Ui 

v' 

=  D 

Vi 

w' 

w' 

(4.15) 


Note  that  in  this  subsection  the  primed  notation  was  only  used  to  be  consistent  with  the 
figure  and  Wilczak  et  al.  (2001).  To  avoid  confusion  with  the  fluctuating  parts,  the  wind 
vector  components  are  redefined  as: 
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(4.16) 


" 

" 

^rot 

u' 

^rot 

= 

v' 

^rot 

w' 

4.4.4  Data  Detrending 

However  the  rotation  procedure  is  not  sufficient  to  prepare  the  data  for  an  analysis 
using  the  auto  covariance  function  since  it  does  not  remove  the  mean  value  from  the  Urot- 
component.  It  also  ignores  any  possible  trend  in  all  three  wind  vector  quantities  Urou  Vrot 
and  Wrot-  Therefore  a  detrending  procedure  has  to  be  applied.  This  is  done  by  fitting  a 
linear  function  through  the  data: 

yfit  =  ax  +  b,  (4.17) 

where  a  is  the  slope  of  the  function  fitted  through  the  data  and  b  is  the  y-intercept.  Here,  x 
is  a  variable,  such  as  time,  and  yja  is  a  variable  defending  on  x,  which  can  be  any  of  the 
velocity  components  or  the  temperature.  The  constants  a  and  b  can  be  determined  as 


a  = 


(4.18) 


and 


b  = 


Ef=i^n  Ef=ii/0  -  Ef=i^0  Ef=i^.i/^ 


(4.19) 


where  N  is  the  number  of  points  of  the  original  data,  i  is  the  index  of  a  point  in  the  original 
data,  Xi  is  the  x-value  of  a  point  with  the  index  i  in  the  original  data  and  yi  is  the  y-value 
of  a  point  with  the  index  i  in  the  original  data.  The  original  data  can  be  described  as 


Vi  = 


(4.20) 
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with 


*  =  1,2,3,  ...iV,  (4.21) 

where  /()  means  function  of.  The  function  yfu  has  now  to  be  subtracted  from  the  original 
function  yielding  a  detrended  function: 


Vdeti  =  f{xi)  -  [axi  +  b]. 


(4.22) 


4.4.5  Example  for  Raw,  Rotated  and  Detrended  Time  Series 

To  illustrate  the  effects  of  rotation  and  detrending  on  the  wind  vector  components  u,  v 
and  w,  a  35  min  time  series  is  shown  as  an  example  in  figure  4.5.  In  figure  4.5  the  diagrams 
in  the  left  column  show  the  raw  time  series.  The  center  column  displays  the  rotated  time 
series  and  the  right  column  contain  first  rotated  and  then  detrended  time  series.  The  time 
series  from  the  top  to  bottom  are  u,  v  and  w,  respectively.  Looking  at  the  diagrams  on  the 
left  hand  side,  it  can  be  seen  that  u  and  v  are  varying  more  strongly  than  w.  The  center 
diagrams  show  that  the  rotation  concept  can  also  be  seen  as  transferring  the  mean  parts 
of  V  and  w  to  the  u  component.  The  over-plotted  lines  in  these  diagrams  are  the  linear 
regression  lines  used  to  detrend  data.  It  can  be  seen  that  there  is  a  slight  upward  trend, 
even  if  the  mean  is  zero  as  shown  for  v  and  w.  In  the  three  diagrams  on  the  right  only 
the  fluctuating  parts  are  left,  since  the  rotation  and  then  detrending  was  applied,  so  that  the 
linear  trend  is  eliminated  and  the  mean  value  is  zero. 

4.5  Data  Analysis 

This  section  will  analyze  and  describe  the  data  from  several  perspectives.  First,  the 
noise  level  and  the  measurement  errors  are  estimated,  in  subsections  4.5.1  and  4.5.2  re¬ 
spectively,  to  provide  an  understanding  of  the  precision  of  the  data.  The  second  step  is  the 
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Raw  w  [m/s]  Raw  v  [m/s]  Raw  u  [m/s] 


2005-07-20  1 5:41 :1 5  -  2005-07-20  16:16:15 
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Figure  4.5.  Wind  vector  components  for  a  35  min  time  series  block  (left  column  -  raw), 
(center  column  -  rotated),  (right  column  -  detrended) 


47 


time-domain  analysis  using  the  autocovariance  function  and  the  structure  function  in  sub¬ 
section  4.5.3.  Subsection  4.5.4  presents  a  spectral  viewpoint  by  means  of  the  periodogram 
and  the  power  spectrum.  Based  on  the  time  and  spectral  domain  analysis,  subsection  4.5.5 
shows  the  inertial  subrange  and  presents  the  velocity  and  temperature  structure  parameter. 
In  subsection  4.5.6  the  integral  length  and  time  scales  are  resolved  and  discussed.  The  final 
subsection  4.5.7  shows  the  approach  of  isotropic  behavior  based  on  structure  function  and 
power  spectrum. 

4.5.1  Noise  Estimation 

A  measured  time  series  m{t)  can  be  described  as  a  superposition  of  a  true  time  series 
f{t)  and  a  noise  component  n{t): 


m{t)  =  fit) +n{t).  (4.23) 

Assuming  that  the  noise  n{t)  is  zero-mean  and  uncorrelated,  its  ACF  becomes: 

[  :  r  =  0 

Bn{r)  =  I  (4.24) 

I  0  :  r^O 

Here  af  is  the  noise  variance.  Now  considering  the  ACF  of  m(f): 

Bmij)  =  {[f{t)+n{t)][f{t+T)+n{t+T)])  =  {f{t)  +  f{t+T))  +  {n{t)+n{t+T)),  (4.25) 

where  statistical  independence  is  assumed,  such  that: 

{f{ti)n{t2))  =  0  (4.26) 
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for  any  ti  and  t2.  Therefore 


Bmir)  =  B  fir)  +  (4.27) 

Now,  considering  the  difference  between  5^(0)  and  Bmir): 

Bm{0)  -  Bmir)  =  Bf{0)  +  -  [5/(t)  +  0]  =  5/(0)  -  Bf{T)  +  (4.28) 

or 

ABUT)  =  ABf{T)  +  al  (4.29) 

where 

ABm{r)  =  Bm{0)  -  Bm{r),  (4.30) 

and 

A5/(r)  =  5/(0)-5/(r).  (4.31) 

If  r  <<  Ti,  then  A5/(r)  may  or  may  not  be  negligible  compared  to  a^.  If  A5/(r) 
negligible,  then  equation  [4.29]  reduces  to 

al^AB^r).  (4.32) 

In  many  cases  A5/(r),  however,  may  not  be  neglected.  Using  equation  [2.18]  A5/(r) 
can  then  be  written  in  terms  of  the  structure  function: 

AB,(t)  =  (4.33) 

If  r  is  in  the  inertial  subrange,  equation  [2.40]  can  be  inserted  into  equation  [4.33],  which 
gives  the  following  result: 

A5/(r)  =  (4.34) 
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Inserting  equation  [4.34]  into  equation  [4.29]  gives: 


=  A5^(r)  -  (4.35) 

valid  only  for  the  inertial  subrange.  For  sonic  anemometer,  an  effect  known  as  spatial 
averaging  smooths  out  all  variations,  which  are  smaller  then  the  measurement  path  d.  Cor¬ 
respondingly,  a  range  of  time  lags  can  be  determined  using  Taylors  hypothesis  in  which  the 
ACF  will  be  low-pass  filtered: 


<i  =  ^  =  0,082s,  (4,36) 

s 

where  d  =  0.15  m  was  found  from  R.M.  Young  Company  (1999)  and  U  =  1.84^  was  com¬ 
puted  from  the  time  series  of  the  velocity  component  u  from  data  set  2005-10-18  16:29:52 
EDT  to  2005-10-18  17:29:52  EDT.  Eigure  4.6  shows  the  ACE’s  of  the  velocity  components 
u,v,w  and  the  temperature  9  in  the  vicinity  of  the  zero  lag  for  the  corresponding  time.  Only 
the  ACE  of  the  temperature  shows  a  spike  at  the  zeroth  lag.  Eooking  to  the  correspond¬ 
ing  detrended  temperature  time  series  in  figure  4.7  shows  that  the  temperature  fluctuation 
during  this  measurement  period  were  in  the  same  range  as  the  precision  of  ±0.1  K  for  the 
temperature  measurements. 

Eor  the  velocity  data,  however,  during  the  whole  experiment  the  fluctuation  were  many 
times  larger  then  the  precision  of  about  ±0.05^  given  by  R.M.  Young  Company  (1999) 
for  velocity,  so  that  a  noise  spike  could  not  be  found.  Eigure  4.8  shows  these  velocity 
fluctuation  correspondingly  to  the  temperature  fluctuation  of  figure  4.7.  The  noise  in  the 
temperature  data  can  be  estimated  by  using  equation  [4.32]  as 
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{T)[{m/sf]  B^{T)[{m/sY\ 


2005-10-18  16:29:52-2005-10-18  17:29:52 


Time  lag/Seconds  Time  lag/Seconds 


Time  lag/Seconds 


Time  lag/Seconds 


Figure  4.6.  ACF  of  u,v,w  and  6  for  maximum  time  lag  r  =  0.32s 
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Figure  4.7.  Time  series  of  detrended  temperature  9 
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10 


2005-1 0-1 8  1 6:29:52  -  2005-1 0-1 8  1 7:29:52 

10| - ^ - 


10 


Figure  4.8.  Time  series  of  rotated  and  detrended  wind  veloeity  components  u,v  and  w 
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alg  ^  AB^g{Ts) 

~  -Bm(?(0)  Buigi^Tg) 

=  O.OOSTir^  -  O.OOSOir^ 

=  O.OOOTir^ 

(4.37) 

or 

an0  =  0.026i^,  (4.38) 

where  is  the  variance  of  the  noise  added  to  the  temperature  time  series,  ane  the  cor¬ 
responding  standard  deviation,  B^g  the  measured  ACF  of  the  temperature  plus  noise  and 
Tg  =  0.0314s  the  sampling  interval,  or  first  time  lag.  Trying  to  use  equation  [4.35]  instead, 
with  U  =  1.84—,  Cn  =  0.02-^,  Tg  =  0.0314s  and  ABmeirs)  =  0.00077f^,  results  in: 

®  m3 

a^g  =  AB^g{Tg)  - 

=  0.000777^  -  ^0.02-^(1.84— 0.0314s)i 
2  77^3  S 

=  0.000777^  -  0.0015077^  =  -0.000877^, 

(4.39) 

which  cannot  be  true  since  there  is  no  negative  noise,  and  it  is  likely  due  to  the  explained 
spatial  averaging  effect. 

4.5.2  Measurement  Errors 

To  describe  measurement  performance  the  specification  of  sensors  usually  contains  the 
terms  accuracy,  precision,  sensitivity  and/or  resolution,  which  are  reviewed  in  the  follow¬ 
ing.  The  precision  P  (or  sensitivity)  is  defined  as  the  standard  error  of  the  fluctuation: 

P  =  -t'Y),  (4.40) 
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where  m'  are  the  measured  fluctuations  and  t'  are  the  true  fluctuations.  Now  consider  the 


standard  error  A  (referred  to  as  accuracy)  of  a  non-zero-mean  variable: 

A  =  \/ ((m  -  t)2) 

=  A/{[((m)  -f  m')  -  {{t)  + 


(4.41) 


The  accuracies  for  the  wind  speed  U  and  the  temperature  6  are  given  in  table  4.3  and 
4.4.  In  subsection  4.5.1  the  noise  in  the  temperature  data  was  estimated  as  ane  =  0.026^^. 
Assuming  a  normal  distribution,  99.7%  of  all  values  are  within  ±3ane  =  ±0.0787^  of  the 
mean  value.  This  can  be  regarded  as  the  precision  of  the  temperature  measurment,  which  is 
much  better  then  the  accuracy  shown  in  table  4.4.  It  can  be  concluded  that  the  mean  value 
of  the  temperature  is  estimated  much  worse  then  the  fluctuations  of  the  temperature.  To 
confirm  the  accuracy  for  temperature  or  velocity  an  external  calibration  standard  would  be 
needed. 

For  the  wind  speed  U  table  4.3  gives  an  equation  to  calculate  the  accuracy  of  the  standard 
deviation  au-  Choosing  =  1-ly  from  the  the  dataset  2005-10-19  23:29:52  EDT  to 
2005-10-20  00:29:52  EDT  gives: 

A  = 


(4.42) 

While  it  was  not  possible  to  estimate  the  noise  in  the  wind  speed  data  in  subsection  4.5.1, 
an  order  of  magnitude  estimate  is  possible  from  the  power  spectmm,  as  performed  later  in 
subsection  4.5.4.  The  standard  deviation  of  the  noise  derived  from  the  power  spectrum  is 


1  m 

± - +  ±0.05— 

100  s 

1  m  m 

± - 1.1—  ±±0.05  — 

100  s  s 

m 

±0.06— 
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0.083y,  which  results  in  a  precision  of  ±0.24^.  This  of  course  is  larger  than  0.06^.  As 
pointed  out  in  subsection  4.5.4  as  well,  this  value  is  overestimated  and  nevertheless  agrees 
by  a  factor  of  about  four  to  five  with  the  range  of  the  accuracy.  A  typical  sensor  contains  an 


Range 

0  to  40  m/s 

Resolution 

0.01  m/s 

Threshold 

0.01  m/s 

Accuracy 

±1  %  rms  ±  0.05  m/s  (0  to  30  m/s) 
±3  %  rms  (30  -  40  m/s) 

Table  4.3.  Performance  on  wind  speed  from  R.M.  Young  Company  (1999) 


Range 

-50° C  to  50°C' 

Resolution 

0.01°^ 

Accuracy 

±2°^  (0  to  30  m/s) 

Table  4.4.  Performance  on  temperature  from  R.M.  Young  Company  (1999) 


analog  to  digital  converter  (ADC)  with  a  certain  number  of  bits  Nadc-  Then  the  resolution 
R  or  the  quantization  step  Q  (assuming  linear  quantization)  is  defined  as: 


R  =  Q 


range 


2N 


(4.43) 


where  range  is  the  allowed  input  range  of  the  sampled  analog  signal,  which  can  be  uni- 
or  bipolar.  Thus  the  value  of  the  resolution  is  not  useful  for  any  performance  description. 
However,  it  should  be  better  than  the  precision  and  the  accuracy,  otherwise  it  will  limit  the 
performance.  It  can  be  seen  from  table  4.3  and  4.4  that  these  requirements  are  fulfilled  for 
the  used  sonic. 
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4.5.3  Autocovariance  Function  and  Structure  Function 


In  figure  4.9  the  ACF  of  the  velocity  component  u  is  shown  for  the  measurement  heights 
of  3.0  m  and  3.8  m.  The  maximum  lag  used  is  r  =  600s.  It  can  be  seen  that  there  are  two 
distinct  time  lag  ranges.  The  first  is  around  the  zero  lag  up  to  about  |r  |  ~  50  —  60s,  and  the 
other  one  is  the  range  beyond  this. 

Following  Orlanski  (1975),  the  time  and  length  scales  of  atmospheric  processes  can  be 
divided  in  categories  of  micro-,  meso-  and  macro-scales.  These  categories  can  further  be 
subdivided  into  subcategories  named  gamma,  beta  and  alpha.  From  figure  4.10,  it  can  be 
seen  that  the  range  in  the  vicinity  of  the  zero  lag  is  denoted  the  micro-gamma  region  with 
time-scales  until  1  min  and  length-scales  until  20  m.  An  extract  of  atmospheric  processes 
classified  by  Orlanski  (1975)  for  the  micro  range  is  shown  in  table  4.5. 

To  extract  the  turbulent  structure  parameter  C|,  the  stmcture  function  Dfij)  is  computed 
assuming  stationarity  from  the  ACF  Bfij)  using  equation  [2.18].  Figure  4. 11  shows  Du{t) 
computed  from  Bu{t)  shown  in  figure  4.9.  Ideally  the  structure  function  approaches,  for 
large  time  lags  r,  a  constant  value  of  2i?j(r)  =  2a‘j.  In  4.11,  the  asymptotic  behavior  is 
approximately  reached  at  about  r  =  Imin. 


Range 

Time  Scale  r 

Length  Scale  R 

Processes 

Micro-Gamma 

r  <  1  min 

R  <  20m 

Plumes 

Roughness 

Turbulence 

Micro-Beta 

Imin  <  T  <5  min 

20m<  R  <  200m 

Dust  Devils 
Thermals 

Wakes 

Micro-Alpha 

5min  <  r  <  30  min 

200m  <  R  <  2000m 

Tornadoes 

Deep  Convection 
Short  Gravity  Waves 

Table  4.5.  Atmospheric  processes  in  the  micro  scale  after  Orlanski  (1975) 
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2005-10-19  23:29:52  -  2005-10-20  00:29:52 


Figure  4.9.  ACF  of  u  with  a  maximum  time  lag  r  =  600s  for  sonic  1  and  2 


58 


200 

m 

20 

Mikro- 

Alpha 

Mikro- 

Beta 

m 

Mikro- 
Gamma  1 

1  1  1  1  1  1  1 

1  min  5  min  30  min  3  h  1  day  3  days 


Figure  4.10.  Scales  of  atmospheric  processes;  after  Orlanski  (1975);  diagram  from  Foken 
(2003)[p.5] 
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2005-1 0-1 9  23 :29 :52  -  2005-1 0-20  0  0  29 :52 


Figure  4.11.  Structure  function  of  u  with  a  maximum  time  lag  r  =  600s  for  sonic  1  and  2 
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4.5.4  Periodogram  and  Power  spectrum 


While  subsection  4.5.3  described  turbulence  from  a  time-domain  perspective,  a  spectral 
description  in  the  form  of  the  power  spectrum  or  the  wavenumber  spectmm  is  often  more 
useful.  Several  ways  are  available  to  estimate  the  one-sided  power  spectrum  S /(cu),  one  of 
which  is  to  use  the  Wiener- Khintchine  relationship  [2.21]  and  equation  [2.43]  based  on  the 
Fourier  transformation  of  the  ACF: 


1  /■“ 

Sfiu)  =  2Wf{u;)  =  2—  /  dr.  (4.44) 

27r 


Another  method  is  the  direct  Fourier  transformation  (FT)  of  the  time  series  to  obtain 
the  spectrum  Lf{u)  as  the  properly  normalized  magnitude  square  of  the  FT: 

Sf{uj)  =  2Wf{uj)  =  2Lf{uj)Lf{u;y  =  2\Lf{uj)f.  (4.45) 


Note  that  the  spectrum  computed  by  this  method  cannot  have  any  imaginary  part,  because 
of  the  absolute  operation  involved. 

The  one-sided  power  spectmm  Sf{uj)  is  also  known  as  the  variance  spectmm,  because  of 
the  following  property: 


(Tj  =  /  Sf{uj)duj,  (4.46) 

Jo 

which  is  the  requirement  for  a  computed  one-sided  power  spectrum  to  be  correctly  nor¬ 
malized.  Table  4.6  shows  such  a  computation  done  for  the  power  spectra  of  u,v,w  and 
9  computed  from  the  dataset  using  both  described  spectral  methods  and  directly  from  the 
time  series.  The  values  shown  agree  well  with  each  other,  as  expected.  The  power  spec¬ 
trum  S'/ (a; )  of  a  random  time  function  should  strictly  be  distinguished  from  the  raw  power 
spectrum  or  the  periodogram,  denoted  Sf{uj).  The  latter  one  is  a  random  function  of  u, 
while  the  power  spectmm  is  a  deterministic  function  of  u. 

However  the  periodogram  can  be  smoothed  to  give  an  estimate  of  the  power  spectrum.  The 
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Variable 

Time  series 

FFT  direct 

FFT  ACF 

1.1971 

1.1971 

1.1970 

0.7611 

0.7611 

0.7610 

BMOil 

0.3169 

0.3169 

0.3168 

MBiB 

0.0407 

0.0407 

0.0407 

Table  4.6.  Variance  computed  from  time  series,  direct  FFT  and  ACF  FFT 


spectral  smoothing  technique  applied  is  logarithmic-equidistant-averaging,  which  averages 
all  data  values  in  a  logarithmic  constant  intervals  of  a  given  width,  as  used  in  Siebert  and 
Muschinski  (2001).  Figure  4.12  shows  an  example  of  the  real  and  imaginary  parts  of  the 
periodogram  and  the  power  spectrum  of  u  computed  using  the  ACF.  Note  that  the  imagi¬ 
nary  part  is  nearly  flat  over  the  spectmm  and  about  10  orders  of  magnitudes  smaller  then 
the  real  part,  which  designates  it  as  numerical  noise  due  to  unavoidable  rounding  errors. 
There  should  indeed  be  no  imaginary  part  due  to  the  computation,  since  the  ACF  is  an  even 
and  real  function,  and  the  Fourier  transformation  of  an  even  and  real  function  is  purely  real. 
Thus  the  results  of  the  real  operation  and  the  absolute  operation  become  identical.  The  pe¬ 
riodogram  and  the  power  spectrum  are  shown  in  figure  4.13  for  the  method  using  the  ACF 
to  compute  the  FFT  and  in  figure  4.14  for  the  method  using  the  time  series  directly.  Now 
comparing  figure  4.13  and  figure  4.14  the  power  spectrum  with  the  /  “§  law  shows  a  good 
agreement  for  nearly  two  decades  from  about  0.2  Hz  to  10  Hz.  However  at  frequencies 
between  10  Hz  and  the  Nyquist  frequency  of  16  Hz  a  slight  flattening  of  the  spectrum  can 
be  seen.  Reasons  for  this  can  be  that  the  noise  level  of  the  sonic  is  reached  or  that  energy 
is  aliasing  back  from  above  the  Nyquist  frequency. 

It  is  possible  to  conservatively  estimate  the  variance  of  the  noise  as  following: 

p^nyq  P^nyq 

al=  Sn{oo)duj  =  Sn  doj,  (4.47) 

Jo  Jo 

where  Unyq  is  the  angular  Nyquist  frequency.  Assuming  a  flat  noise  floor  and  that  the 
smoothed  measured  power  spectmm  is  severely  limited  by  noise  at  the  Nyquist  frequency 


62 


2005-1 0-1 9  23:29:52  -  2005-1 0-20  00:29:52 


Freq/Hz 


Figure  4.12.  Real  and  imaginary  parts  of  periodogram  and  power  speetrum  of  u  time  series 
eomputed  from  ACF 
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2005-10-19  23:29:52  -  2005-10-20  00:29:52 
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Figure  4.13.  Periodogram  and  power  speetrum  of  u  time  series  eomputed  from  ACF 
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vjrs 


2005-10-19  23:29:52  -  2005-10-20  00:29:52 


Periodogram 

Estimated  Powerspectrum 


Figure  4.14.  Periodogram  and  power  spectrum  of  u  time  series  computed  directly  from 
time  series 
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equation  [4.47]  can  be  written  as: 


(Jn  —  S  f{uJnyq)^^- 


nyq- 


(4.48) 


Table  4.7  lists  the  variance  of  the  noise  for  all  velocity  components  and  for  temperature 


Variable 

FFT  direct 

FFT  ACF 

ACF 

0.0069 

0.0069 

N/A 

0.0088 

0.0088 

N/A 

BMIHil 

0.0102 

0.0103 

N/A 

— BM 

0.00069 

0.00069 

0.0007 

Table  4.7.  Variance  of  noise  computed  from  power  spectrum  and  ACF 


Variable 

FFT  direct 

FFT  ACF 

ACF 

Siebert  and  Muschinski  (2001) 

1 

0.083 

0.083 

N/A 

0.02 

S 

1 

0.093 

0.093 

N/A 

0.02 

0.101 

0.101 

N/A 

0.02 

CTnO  [K] 

0.026 

0.026 

0.026 

0.02 

Table  4.8.  Standard  deviation  of  noise  computed  from  power  spectrum  and  ACF 


using  the  spectral  techniques.  Table  4.8  shows  the  corresponding  standard  deviations  and, 
in  addition,  values  measured  by  Siebert  and  Muschinski  (2001)  with  a  similar  sonic.  It 
can  be  seen  that  the  variance  of  the  noise  in  the  temperature  data  agrees  well  with  all  tree 
methods,  and  is  close  to  the  value  from  Siebert  and  Muschinski  (2001).  The  variance  of 
the  noise  in  the  velocity  components  might  have  been  overestimated  because  the  flattening 
in  the  power  spectra  is  not  very  strong  and  there  is  no  visual  noise  spike  in  the  ACF  of  the 
velocity  components.  Comparing  the  velocity  components  to  the  smaller  values  retrieved 
by  Siebert  and  Muschinski  (2001)  supports  this  interpretation. 

In  figure  4.15  the  power  spectra  of  the  velocity  components  are  computed  with  the  two 
described  spectral  methods  and  plotted.  It  can  be  seen  that  both  spectral  analysis  methods 
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lead  to  very  similar  results.  Further,  it  ean  be  noted  that  all  obey  the  5/3  law  over  a  range 
of  nearly  two  deeades.  However  S'^(/)  has  a  signifieant  lower  power  density  at  lower 
frequeneies,  which  can  be  explained  by  the  presence  of  the  ground.  The  ground  prevents 
large-scale  or  long-term  vertical  motion  at  low  heights  (assuming  the  ground  is  horizontal 
and  flat).  Similarly  in  figure  4.16  the  temperature  power  spectmm  is  plotted  for  both  spec- 


2005-10-19  23:29:52  -  2005-10-20  00:29:52 


Figure  4.15.  Power  spectrum  of  velocity  component  time  series  computed  directly  and 
from  ACF 

tral  methods,  which  again  show  agreement  and  follow  the  5/3  law.  The  effect  of  the  noise 
floor,  however,  is  much  more  visible,  compared  to  the  velocity  power  spectra.  This  is  indi¬ 
cated  by  a  stronger  flattening  of  the  power  spectrum  near  the  nyquist  frequency,  caused  by 
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Seif)[ 


the  lower  signal-to-noise  ratio  (SNR)  of  the  temperature  measurements  compared  to  that 
of  the  velocity. 


2005-10-19  23:29:52  -  2005-10-20  00:29:52 


Figure  4.16.  Power  spectrum  of  temperature  time  series  computed  directly  and  from  ACF 


4.5.5  Inertial  Subrange  and  Structure  Parameter 

In  subsection  4.5.3  the  stmcture  function  Df{T)  was  computed  and  an  example  for 
^n('r)  was  given  in  figure  4.11.  The  turbulent  structure  parameter  Cj  can  be  computed 
from  Df=u{'T')  based  on  equation  [2.40]: 


Df  jr) 

f/2/3^2/3  • 


(4.49) 
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In  figure  4.17  the  velocity  structure  parameters  C'^(r),C^(r)  and  C'.^(r)  are  shown  while  in 
figure  4.18  the  temperature  structure  parameter  Cg  (r)  is  plotted.  It  needs  to  be  kept  in  mind 
that  equation  [4.49]  is  valid  only  for  the  inertial  subrange.  This  range  can  be  recognized 
as  a  plateau.  In  addition,  the  25^(0)  asymptote  is  plotted.  The  intersection  of  the  2/3  law 
and  the  25^(0)  asymptote  can  be  defined  as  the  intersection  time  lag  tq  (here  Tu,Tv,Tu,  and 
T0  respectively).  Values  of  the  structure  parameters  can  also  be  derived  from  the  power 


2005-10-19  23:29:52  -  2005-10-20  00:29:52 


Figure  4.17.  Structure  parameter  of  velocity  components  from  structure  functions 
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=  Dfl(r)/(f7T)2/3 


2005-10-19  23:29:52  -  2005-10-20  00:29:52 


Figure  4.18.  Structure  parameter  of  temperature  Cg  (r)  from  structure  function 
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spectra  by  rewriting  equation  [2.45]  as: 


cjif) 


4.02Sf{f) 


(4.50) 


Figure  4.19  shows  the  results  from  equation  [4.50]  for  the  velocity  structure  parameter  and 
figure  4.20  shows  the  structure  parameter  of  the  temperature  derived  from  the  power  spec¬ 
trum. 

Table  4.9  lists  the  structure  parameter  derived  with  the  spectral  methods  and  the  structure 
function.  It  can  be  seen  that  in  all  cases  is  the  largest,  followed  by  Cl  and  then 
While  both  spectral  methods  agree  closely  between  themselves,  the  values  extracted  from 
the  stmcture  function  deviate  from  them.  In  the  case  of  deriving  the  stmcture  parameter 
from  the  structure  function  only  a  peak  instead  of  a  well  defined  plateau  is  indicating  the 
2/3  law,  and  therefore  the  inertial  subrange.  Therefore  structure  parameter  can  only  be 
chosen  as  the  maximum.  However,  this  is  not  the  case  when  it  is  derived  from  the  spectral 
methods,  in  which  case  the  inertial  subrange  is  wider  and  can  be  recognized  as  a  plateau  in 
figure  4.19  and  4.20  at  higher  frequencies. 

Comparing  these  two  ways  (spectral  methods  and  structure  function),  the  way  using  the 
spectral  methods  provides  a  much  wider  inertial  subrange,  which  allows  for  a  better  es¬ 
timate  of  the  structure  parameter,  rather  than  just  using  the  peak  value  from  the  structure 
function  estimate. 


Variable 

FFT  direct 

FFT  ACF 

Structure  Function 

(  {m/sY 
t  rn,2/3 

] 

0.4529 

0.4532 

0.3415 

/  (m/sY 

V  rf7,2/3 

] 

0.5506 

0.5512 

0.3823 

Y  {■m./sY 

A  m?/^ 

] 

0.3868 

0.3869 

0.3145 

\  K  '^ 

0.0173 

0.0173 

0.0154 

Table  4.9.  computed  from  power  spectra  and  structure  function 
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2005-10-19  23:29:52  -  2005-10-20  00:29:52 


Figure  4.19.  Structure  parameter  of  velocity  components  from  power  spectrum 


4.5.6  Integral  Time  and  Length  Scales 

The  intersection  time  lag  tq  given  in  subsection  4.5.5  can  be  determined  from  the  ACF 
by  assuming  a  2/3  dependence  in  following  form: 

B,(t)  =  B,(Q)  -  ]^D,  =  a)-  (4.51) 

where  equation  [2.40]  was  inserted  into  equation  [2.18]  to  resolve  equation  [4.51].  To  find 
the  intersection  point  of  the  25/(0) -asymptote  and  the  2/3  law,  equation  [4.51]  can  be 
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{I)oszo-f  =  (Dio 


2005-10-19  23:29:52  -  2005-10-20  00:29:52 


Figure  4.20.  Structure  parameter  of  temperature  from  power  spectmm 


rewritten  as: 


2a  j  =  CjU^T^, 


(4.52) 


where  the  time  lag  r  was  set  to  the  intersecting  time  lag  tq.  The  intersection  time  lag  tq 
can  now  be  written  as: 


To  =  25 


[cpu 


(4.53) 


Using  equation  [2.31],  insert  [4.51]  and  assuming  i?/(r)  =  0  for  r  >  r*,  r*  can  be  written 


as: 
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Ti  = 


5/(0)  Jo 


=  2r, 


=  2 


5/(0)  Jo 

CjU^/^3  S 


5/(r) dr 


Tn 


a. 


3 

0  0  d 

3 

3  o-j 

0  0  4 

[C]]lu_ 

aj  5 

[Cj]lu_ 

5  3  5 

22  -  -22 


O'/ 


5  J  [CpU 


=  r^O- 


(4.54) 


Figures  4.21,  4.22,  4.23  and  4.24  show  the  application  of  this  technique,  from  here  it  can  be 
seen  that  the  analytical  result  for  the  relationship  r*  =  |ro  from  equation  [4.54]  is  confirmed 
by  the  numerical  integration  over  (r).  Using  Taylor’s  hypothesis  the  integral  time  scale 
Ti  can  be  converted  into  the  integral  length  scale  rf. 

ri  =  UTi  (4.55) 


and  correspondingly  for  the  integral  scales  of  the  angular  frequency  oji  (spectral  width)  and 
the  wavenumber  k*: 

oJi  =  Uki.  (4.56) 

Table  4.10  shows  the  results  of  the  integral  scales  calculated  using  equations  [4.55],  [4.56], 
[2.33]  and  [2.36].  Two  effects  on  the  length  scales  can  be  seen  in  table  4.10.  First  the  length 
scale  Tiy^  =  1.93  m  of  tu  is  much  smaller  than  the  length  scale  =  12.10  m  of  m  and  the 
length  scale  =  6.00  m  of  v.  Secondly  the  length  scale  of  w  is  about  1/3  less  than  the 
height  =  3m,  while  the  length  scales  of  u  and  v  are  twice  and  four  times  the  height  2: 
respectively.  Looking  at  ACF’s  shown  in  figure  4.23  for  w,  in  figure  4.21  for  u  and  4.22 
for  V,  a  faster  decline  of  the  ACF  for  w  can  be  seen  then  for  u  and  v.  Another  observation 
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Var. 

Time  Scale 

Length  Scale 

Spectral  Width 

Spectral  Width 

Spectral  Width 

n 

Ti 

LOi 

h 

Ki 

[S] 

[m] 

[^] 

[Hz] 

[-] 

u 

4.99 

12.10 

1.26 

0.20 

0.52 

V 

2.48 

6.00 

2.53 

0.40 

1.05 

w 

0.80 

1.93 

7.76 

1.23 

3.26 

e 

3.81 

9.24 

1.64 

0.26 

0.68 

Table  4.10.  Integral  scales 


is  that  u  has  the  largest  length  scale.  This  is  consistent  with  the  idea  that  the  ground  limits 
the  length  scale  of  w  as  mentioned  before. 

Another  way  to  analyze  the  length  scale  and  the  time  scale  is  to  compute  their  values  for 
l-hour  blocks  over  the  whole  experiment  time  and  to  compare  their  behavior  to  the  mean 
wind  speed  computed  over  the  same  time  interval.  Figure  4.25  shows  the  mean  wind  speed 
U  and  the  integral  time  scale  Tiu.  Using  equation  [2.39]  the  integral  length  scale  riu  as  a 
function  of  time  is  plotted  in  figure  4.26.  Calculating  the  mean  value  fi^  and  the  median 
value  median{riu)  of  the  integral  length  scale  r^u  shown  in  figure  4.26  results  in: 

fiu  =  12.05m  (4.57) 


and 


medianir'iu)  =  10.67m. 


(4.58) 


It  can  be  observed  that  U  and  Tiu  are  inversely  correlated  to  each  other  in  figure  4.25,  while 
in  figure  4.26  it  can  be  seen  that  the  integral  length  scale  has  a  much  more  constant  behavior. 
From  this,  it  can  be  said  that  the  size  of  the  stream  wise  mean  velocity  length  scale  does 
not  dependent  much  on  the  wind  speed.  Instead  it  can  be  found  in  Muschinski  and  Roth 
(1993)  and  Muschinski  et  al.  (2004)  that  it  depends  on  the  height  and  it  is  known  as  the 
largest  measurable  wavelength  or  the  isotropic  cut-off  wavelength  Xiu.  It  was  found  to  be: 


A,',,  =  4z. 


(4.59) 
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Using  a  height  z=3m  results  in  a  wavelength  of  12ni,  whieh  is  a  elose  mateh  to  the  mean 
value  and  of  the  same  order  as  the  median  value  eomputed  before.  Notice,  that  strictly 
speaking  the  integral  length  scale  riu  and  the  isotropic  cut-off  wavelength  are  not  nec¬ 
essarily  identical,  but  they  should  be  of  the  same  order  of  magnitude. 


2005-10-19  23:29:52  -  2005-10-20  00:29:52 


Figure  4.21.  Autocovariance  function  Bu{t)  and  integral  time  scale 
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2005-10-19  23:29:52  -  2005-10-20  00:29:52 
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Figure  4.22.  Autocovariance  function  B^{t)  and  integral  time  seale 
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2005-10-19  23:29:52  -  2005-10-20  00:29:52 


Figure  4.23.  Autocovariance  function  and  integral  time  seale  Tiy^ 
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2005-10-19  23:29:52  -  2005-10-20  00:29:52 
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Figure  4.24.  Autocovariance  function  Beij)  and  integral  time  scale  Tie 
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U{t)/{m/ s),Tiu{t)/ 


2005-1 0-1 7  1 7:59:52  -  2005-1 0-20  1 1 :59:52 


Figure  4.25.  Integral  time  seale  Tiu  and  mean  wind  speed  U 
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2005-1 0-1 7  1 7:59:52  -  2005-1 0-20  1 1 :59:52 


Figure  4.26.  Integral  length  scale  Viu 
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4.5.7  Isotropy 


Kolmogorov  (1941a)  predicts  a  4/3  ratio  between  the  longitudinal  and  the  transverse 
correlation  functions  of  the  velocity  in  the  inertial  subrange  for  homogeneous  and  isotropic 
turbulence.  Using  equation  [2.54]  and  Taylor’s  hypothesis,  equation  [2.39]  can  be  written 


as: 


and 


Therefore, 


Dw{t)  _  Dv{t)  _  4 
Du{t)  Du{t)  3’ 


Dwjr) 

D.ir) 


Inserting  equation  [4.49]  and  equation  [4.50]  accordingly: 


(4.60) 


(4.61) 


(4.62) 


D^{t)  Sjj{u) 


and 

Dw{t)  ^  Dv{t)  ^  Sy,{uj)  ^  Syjuj) 

Du{t)  Du{t)  Su{uj)  Su{u})' 

Figure  4.27  shows  the  ratios  of  the  structure  functions.  In  figure  4.28  the  ratios  of  the 
power  spectra  are  shown. 

In  figure  4.27  it  can  be  seen  that  the  4/3  ratio  is  not  exactly  reached  for  the  smallest  time 
lag,  which  could  be  a  result  of  spatial  filtering  (as  a  result  of  the  finite  size  of  the  sonic)  that 
attenuates  the  small  scale  turbulence.  The  sharper  decline  of  compared  to  is  a 
likely  result  of  anisotropy  of  the  w  velocity  component,  due  to  the  presence  of  the  ground. 
The  velocity  components  u  and  v,  on  the  other  hand,  have  more  space  available  than  the 
velocity  component  w,  since  they  are  oriented  into  the  horizontal  plane.  Thus  the  isotropic 
behavior  shown  by  is  more  extended  than  the  one  of  444  •  A  similar  behavior  can 

•'  Du(t)  Du{t) 

S  ( 

be  seen  in  figure  4.28.  The  spike  in  figure  4.28  in  the  graph  of  at  /  ^  0.02  Hz  and 
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Figure  4.27.  Ratios  between  strueture  funetions 
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Sn,{f)  S^if)  54/) 


2005-10-19  23:29:52  -  2005-10-20  00:29:52 


Figure  4.28.  Ratios  between  power  spectra 
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in  figure  4.27  in  graph  at  r  1.1  s  seems  to  be  from  a  larger  scale  effect  that  is  not 
directly  related  to  small-scale  turbulence. 
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CHAPTER  5 


SUMMARY  AND  CONCLUSIONS 


5.1  Summary 

In  chapter  one,  motivation  for  turbulence  measurements  in  the  surface  layer  was  given 
and  it  was  pointed  out  that  the  sonic  anemometer  is  a  good  tool  for  such  measurements. 
The  historical  development  of  turbulence  research  was  outlined  and  the  current  state  sum¬ 
marized.  Chapter  two  reviewed  the  theory  of  locally  homogeneous  and  isotropic  turbulence 
and  the  Monin- Obukhov  similarity  theory.  In  chapter  three  the  measurement  system  was 
described.  In  chapter  four,  66  hours  of  wind  velocity  and  temperature  data  including  two 
diurnal  cycles  retrieved  by  means  of  two  sonic  anemometers  and  thermometers  placed  in 
the  atmospheric  surface  layer  at  heights  of  3  m  and  3.8  m  have  been  analyzed.  The  data 
was  divided  into  1  hour  blocks.  Following  a  rotation  from  the  sonic  coordinate  system  to 
a  stream-wise  coordinate  system,  a  detrending  algorithm  was  applied  to  the  data,  so  that 
only  the  fluctuating  parts  were  left.  The  periodograms  of  the  wind  velocity  components, 
Suif),  Sy{f),  and  the  periodogram  of  the  temperature,  ^e(/),  were  resolved  and 

presented  by  applying  the  FFT  to  the  time  series  directly  on  the  one  hand  and  by  comput¬ 
ing  the  autocovariance  function  and  using  the  Wiener-Khintchine  relationship  on  the  other 
hand.  By  means  of  logarithmic-equidistant-averaging  of  the  periodograms  the  power  spec¬ 
tra  of  the  wind  velocity  components,  5'„(/),  Sy{f),  S'^(/),  and  the  power  spectrum  of  the 
temperature,  S0{f),  were  estimated.  The  noise  variances  of  the  wind  velocity  components, 
aln,  ^wn  Were  estimated  from  the  velocity  spectra.  The  noise  variance  of  the  temper- 
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ature,  was  estimated  from  the  power  speetra  and  the  autoeovarianee  funetion.  These 
noise  varianees  were  then  used  to  eompute  the  preeision  of  the  veloeity  and  the  temperature 
measurements.  The  temporal  strueture  funetions  T*«(r),  Dw{t),  D0{t)  and  spatial 

strueture  funetions,  T>«(r),  Dy{r),  D^{r),  Dg{r)  of  the  veloeity  and  temperature  were  eom- 
puted.  These  and  the  power  speetra  were  used  (as  two  different  ways)  to  estimate  the  stme- 
ture  parameters  of  the  wind  veloeity  eomponents  and  the  strueture  parameter 

Cg  of  the  temperature.  The  integral  time  seales  Tiu,  Tiy,  of  the  wind  veloeity  eomponents 
and  the  integral  time  seale  Tig  of  the  temperature  were  eomputed  from  the  autoeovarianee 
funetion.  From  the  integral  time  seales,  the  integral  length  seales  and  the  speetral  widths 
were  eomputed.  The  isotropie  eut-off  wavelength  Ac  was  estimated  using  the  mean  and 
the  median  of  the  integral  length  seale  riu  derived  for  all  66  hours  of  the  data.  Finally,  the 
ratios  of  the  strueture  funetions,  ratios  of  the  frequeney  speetra 

s^(f)  ’  plotted  against  time  lag  and  frequeney,  respeetively. 

5.2  Conclusions 

We  offer  the  following  eonelusions: 

(1)  All  power  speetra  follow  a  elear  -5/3  law  for  nearly  two  deeades,  from  about  0.2  Hz 
to  the  Nyquist  frequeney  of  about  16  Hz,  whieh  is  eonsistent  with  Kolmogorov-type  small 
seale  turbulenee.  This  finding  is  in  agreement  with  earlier  observations  (e.g.  Kaimal  et  al. 
(1972),  Musehinski  et  al.  (2004)).  Below  about  0.02  Hz  the  veloeity  speetra  beeome  flat.  It 
ean  be  seen  that  in  this  range  the  speetral  power  density  of  the  vertieal  wind  eomponent  w  is 
about  10  times  lower  than  the  speetral  power  densities  of  the  horizontal  wind  eomponents 
u  and  V,  whieh  ean  be  explained  by  the  suppression  of  w  fluetuations  due  to  the  presenee 
of  the  ground. 

(2)  It  was  found  that  our  preeisions  (here  defined  as  ±3a  to  inelude  99.7%  of  gaussian 
distributed  errors)  eomputed  from  the  power  speetra  for  the  wind  veloeity  eomponents  (u: 
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±0.249^,  v:  ±0.279^,  w:  ±0.303^)  agree  within  an  order  of  magnitude  to  the  speci¬ 
fication  of  the  accuracy  (~  ±0.1^)  given  by  the  manufacturer.  However,  our  precision 
estimation  is  probabily  too  conservative.  For  the  precision  of  the  temperature  measurement 
it  was  found  that  the  estimation  based  on  the  power  spectrum  and  the  autocovariance  func¬ 
tion  agree  very  closely  to  each  other.  It  was  found  to  be  about  ±0.078  K  (about  25  times 
better  then  the  specified  accuracy  of  ±2  K),  which  is  comparable  to  (±0.06  K)  as  found  by 
Siebert  and  Muschinski  (2001)  for  a  similar,  but  different  type  of  sonic. 

(3)  The  integral  length  scales  (rj„  =  12.1  m,  =  6.0  m,  =  1.93  m,  tie  =  9.24  m)  and 
the  integral  time  scales  (ri„  =  4.99  s,  Tiv  =  2.48  s,  =  0.81  s,  Tie  =  3.82  s)  were  found 
to  be  consistent  with  the  scale  classification  given  by  Orlanski  (1975)  for  atmospheric  tur¬ 
bulence,  (r*  <  20  m  and  r*  <  60  s).  It  was  found  that  the  length  scale  r^u  is  largest,  as 
expected.  Furthermore,  it  was  found  that  the  length  scale  is  on  the  order  of  the  height 
2:  =  3  m  reduced  by  the  roughness  length  zq. 

(4)  The  isotropic  cut-off  wavelength  was  estimated  as  12.05  m  using  the  mean  of  the  in¬ 
tegral  length  scale  of  u  and  as  10.67  m  using  the  median  of  the  integral  length  scale  of  u. 
These  are  close  to  the  theoretical  estimate  from  Muschinski  and  Roth  (1993)  and  Muschin¬ 
ski  et  al.  (2004)  given  as  Ac  ~  4^  ~  12  m. 

(5)  The  temperature  and  velocity  structure  parameters  computed  from  the  structure  func¬ 
tions  and  the  power  spectra  were  found  to  agree  with  each  other  within  about  35%.  A  typ¬ 
ical  expected  range  is  given  by  Ishimam  (1978)  on  [p.344]  as  =  10“®^^^  to  10“^^^^ 
for  the  temperature  structure  parameter  and  =  10“^^^^  to  1^^^  for  the  velocity  struc¬ 
ture  parameter.  Similar  orders  of  magnitudes  for  the  structure  parameters  were  obtained  by 
Muschinski  et  al.  (2004). 

(6)  A  ratio  close  to  4/3  was  observed  for  the  ratio  of  the  structure  functions  for  time 

(,7"  ) 

S  ( 

lags  T  <  0.06  s  and  for  the  ratio  of  the  frequency  spectra  for  frequencies  /  >  2  Hz. 
However,  for  the  ratio  of  the  structure  functions  and  for  the  ratio  of  the  frequency 
spectra  the  approach  to  the  4/3  ratio  is  present  only  at  the  smallest  shown  time  lag 
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and  for  the  highest  shown  frequency.  The  reason  for  this  behavior,  is  the  presence  of  the 
ground,  which  affects  w  more  than  u  and  v. 


5.3  Recommendation  and  Future  Work 

For  future  work,  the  sonics  should  be  used  on  more  flat  and  homogeneous  surface,  to  al¬ 
low  for  more  detailed  statistical  investigations.  Some  possible  future  experiments  include 
the  installations  of  several  sonics  placed  at  the  same  height,  for  the  test  of  Taylor’s  hy¬ 
potheses,  by  comparing  the  temporal  and  the  spatial  ACF  or  structure  functions.  Another 
experiment  would  place  several  sonics  at  significant  different  heights  to  further  confirm  the 
isotropic  cut-off  wavelength  and  to  retrieve  height  profiles  of  velocity,  temperature  and  the 
turbulent  fluxes. 
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